Indices:
1 -> Upstream
2 -> At entry to rotor
3 -> After rotor

Notation:
P -> Pressure
A -> Effective cross-section in static frame
B -> Effective cross-section in rotating frame
c -> Gas velocity in static frame
w -> Gas velocity in rotating frame
α -> Gas deviation angle from radial flow in static frame
β -> Gas deviation angle from radial flow in rotating frame
r -> Distance from axis
ω  -> Angular velocity of turbine
h -> Blade height
a -> Velocity of sound

Necessary power: 9.2 kW
Critical temperature for (Subscript[O, 2]): 155 K

Units:
All units are SI.

\(\(\(eqthermo = {
P\_3 -> P\_C, \  (*\ 
            Chamber\ pressure\ at\ outlet\ *) 

 (*\ 
            No\ expansion\ in\ rotor\ stage\ *) 
P\_2 -> P\_3, \ 

 (*\ 
            Isentropic\ flow\ in\ static\ stage, \ 
            assuming\ A\_1\ is\ infinite\ *) 
M\_2 -> \@\(\(2\/\(γ - 1\)\) \((\((P\_1\/P\_2)\)\^\(\(γ - 1\)\/γ\) - 1)\)\), 
T\_2 -> \(T\_1\) \((1 + \(\((γ - 1)\)\/2\) M\_2\^2)\)\^\(-1\), 
A\_2 -> 
            m\& . \/\(\(c\_2\) ρ\_2\), 

 (*\ 
            Velocity\ of\ sound\ *) 
a\_2 -> \@\(γ\ \(R\_g\) T\_2\), 
c\_2 -> \(a\_2\) 
              M\_2, 

 (*\ 
            Ideal\ gas\ law\ *) 
ρ\_2 -> 
            P\_2\/\(\(R\_g\) T\_2\)
};\)\(
\)
    \)\)
\(\(\(eqgeom\  = \ {
 (*\ 
            Velocity\ triangles\ *) 
α\_2 -> 
            ArcSin[\(ω\ r\_2\)\/c\_2], 
α\_3 -> ArcSin[\(\(w\_3\) Sin[β\_3] + ω\ r\_3\)\/c\_3], 
w\_2 -> \(c\_2\) 
              Cos[α\_2], 
c\_3 -> \@\(\((\(w\_3\) Cos[β\_3])\)\^2 + \((\(w\_3\) Sin[β\_3] + ω\ r\_3)\)\^2\), 

 (*\ 
            Work\ extracted\ from\ turbine\ *) 
W -> 
            m\& . \ \((\(c\_2\) \(r\_2\) Sin[α\_2] - \(c\_3\) \(r\_3\) 
                    Sin[α\_3])\) ω, 

 (*\ 
            No\ expansion\ in\ rotor\ *) 
B\_3 -> 
            B\_2, 
w\_3 -> 
            w\_2, 

 (*\ 
            Effective\ areas\ *) 
B\_2 -> 
            2\ π\ r\_2\ h\_2, 
h\_2 -> \ 
            A\_2\/\(2  π\ \(r\_2\) Cos[α\_2]\), 
h\_3 -> 
            B\_3\/\(2  π\ \(r\_3\) Cos[β\_3]\), 

 (*\ 
            Some\ useful\ numbers\ *) 
Δh\_g -> \(γ\/\(γ - 1\)\) \(m\& . \) \(R\_g\) \((T\_1 - 
                  T\_2)\)
};\)\(
\)
    \)\)
\(\(setvars = {
γ -> 
            1.4, 
R\_g -> 
            8.314/32*^-3, 
P\_C -> 
            3.5*^6, 
P\_\(\(1\)\(\ \)\) -> 
            6*^6, 
m\& .  ->  .44, 
ω -> 2  N[π]\ 85000/60\  (*\ 
            85000\ RPM\ in\ rad/s\ *) , 
T\_1 -> 
            260. , 
r\_2 -> 
            3.254*^-2, 

 (*\ 
            Guesses\ for\ now\ *) 
r\_3 -> 
            1*^-2, 
β\_3 -> \(-28.823239045597287\), 
β\_2 -> 
            0
};\)
\)
alleq = Union[eqthermo, eqgeom];

variableHeads = {Subscript, 
          OverDot};

getArgs[x_] := 
      If[AtomQ[x], x, 
        If[MemberQ[variableHeads, Head[x]], x, 
          Flatten[Map[getArgs, 
              Apply[List, 
                x]]]]]

getVars[x_] := 
      Variables[getArgs[x]]

Print["Variables:"]
vars = getVars[alleq]

  
allsolved = 
        Map[# -> (# //. alleq) &, 
          vars];

  
replaceTarget[x_, r_] := 
      x /. (y_ -> 
              z_) :> (y -> (z /. 
                  r))

  
Print["Solution:"]
numsolved = replaceTarget[allsolved, setvars]

\(\(Δh\_g/\((T\_1 - T\_2)\)\)/\(m\& . \)\  /. 
      numsolved\)

\(\(40000/\((\(Δh\_g/\((T\_1 - 
                    T\_2)\)\)/\(m\& . \))\)\)/\(m\& . \)\  /. 
      numsolved\)

(Subscript[A, 2])*1*^6  /.  numsolved

(Subscript[A, 2])/(2*Pi*(Subscript[r, 2]))  /.  numsolved

(Subscript[α, 2]) /Degree /. numsolved

ω (Subscript[r, 2])  /. numsolved

PlotVsBeta[v_] := 
      Plot[Evaluate[v /. numsolved /. (Subscript[β, 3]) -> 
              x Degree], {x, -90, 90}]

PlotVsOmega[v_] := 
      Plot[Evaluate[v /. numsolved /. ω -> 
              2 π x /60], {x, 0, 1.5*85000}]

PlotVsOmega[W]

\(PlotVsBeta[W] 
      PlotVsBeta[\(1\/2\) \(m\& . \) \((c\_2\^2 - c\_3\^2)\) - 
          W + \(1\/2\) \(m\& . \) \(ω\^2\) \((r\_2\^2 - r\_3\^2)\)] 
      PlotVsBeta[α\_3/Degree] 
      PlotVsBeta[{c\_3\ Sin[α\_3], c\_3}]\)

(Subscript[β, 3])/
        Degree  /.  Solve[
          Evaluate[(Subscript[α, 3]) /. numsolved] == 
            0, (Subscript[β, 3])][[2]]

Old stuff

gV[x_[y__] /; ! MemberQ[variableHeads, x]] := 
      Variables[Flatten[Map[gV, {y}]]]
gV[x_ ] := x

(*eqgeom = {
 (* 
              No expansion in rotor stage *) 
(Subscript[B, 2]) == (Subscript[B, 3]),  

 (* 
              velocity triangle *) 
       (*(Subscript[c, 2]) Cos[(Subscript[α, 2])] == (Subscript[w, 2]) 
                  Cos[(Subscript[β, 2])],*) 
(Subscript[c, 2]) 
                Cos[(Subscript[α, 2])] == (Subscript[w, 2]),   (* 
              This assumes (Subscript[β, 2]) = 
                0 *) 
       (*(Subscript[c, 2]) 
                  Sin[(Subscript[α, 2])] == (Subscript[w, 2]) 
                    Sin[(Subscript[β, 2])] + ω (Subscript[r, 2]),*) 
(Subscript[c, 2]) Sin[(Subscript[α, 2])] == ω (Subscript[r, 2])  (* 
              This assumes (Subscript[β, 2]) = 
                0 *) 

 (*(Subscript[c, 3  r]) == (Subscript[w, 3]) 
                  Cos[(Subscript[β, 3])], 
(Subscript[c, 3  θ]) == (Subscript[w, 3]) 
                    Sin[(Subscript[β, 3])] + ω (Subscript[r, 3]), *) 

 (* 
              Geometric relations on area *) 
 (*(Subscript[A, 2]) == ((2  π (Subscript[r, 2]) (Subscript[h, 2]))/Cos)[(Subscript[α, 2])]*) 
 (*(Subscript[B, 2]) ==  ((2  π (Subscript[r, 2]) (Subscript[h, 2]))/Cos)[(Subscript[β, 2])],*) 
 (*(Subscript[B, 2]) ==  2  π (Subscript[r, 2]) (Subscript[h, 2])*)   (* 
              This assumes (Subscript[β, 2]) = 
                0 *) 
 (*(Subscript[B, 3]) == ((2  π (Subscript[r, 3]) (Subscript[h, 3]))/Cos)[(Subscript[β, 3])]*) 
};*)

\(\(eqthermo = {
P\_3 == P\_C, \  (*\ 
            Chamber\ pressure\ at\ outlet\ *) 

 (*\ 
            No\ expansion\ in\ rotor\ stage\ *) 
P\_2 == P\_3, \ 

 (*\ 
            Isentropic\ flow\ in\ static\ stage, \ 
            assuming\ A\_1\ is\ infinite\ *) 
M\_2 == \@\(\(2\/\(γ - 1\)\) \((\((P\_1\/P\_2)\)\^\(\(γ - 1\)\/γ\) - 1)\)\), 
T\_2 == \(T\_1\) \((1 + \(\((γ - 1)\)\/2\) M\_2\^2)\)\^\(-1\), 
\(m\& . \) == \(A\_2\) \(c\_2\) ρ\_2, 

 (*\ 
            Velocity\ of\ sound\ *) 
a\_2 == \@\(γ\ \(R\_g\) T\_2\), 
c\_2 == \(a\_2\) 
              M\_2, 

 (*\ 
            Work\ *) 
 (*\(W == \(m\& . \)\ \((\(c\_2\) \(r\_2\) Sin[α\_2] - \(c\_3\) \(r\_3\) 
                      Sin[α\_3])\) ω\)\(,\)*) 

 (*\ 
            Ideal\ gas\ law\ *) 
P\_2\  == \ \(ρ\_2\) \(R\_g\) 
              T\_2
};\)
\)
\(eqgeomsolved\  = \ {
α\_2 -> 
          ArcSin[\(ω\ r\_2\)\/c\_2], 
w\_2 -> \(c\_2\) Cos[α\_2], 
w\_3 -> \(w\_2
h\_2 -> \ \(\(A\_2\) Cos[α\_2]\)\/\(2  π\ r\_2\)\), 
B\_2 -> 
          2\ π\ r\_2\ h\_2, 
B\_3 -> 
          B\_2, 
W -> \(m\& . \)\ \((\(c\_2\) \(r\_2\) 
                  Sin[α\_2] - \(c\_3\) \(r\_3\) 
                  Sin[α\_3])\) ω, 
c\_3 -> \@\(\((\(w\_3\) Cos[β\_3])\)\^2 + \((\(w\_3\) Sin[β\_3] + ω\ r\_3)\)\^2\), 

α\_3 -> 
          
          ArcSin[\(\(w\_3\) Sin[β\_3] + ω\ r\_3\)\/c\_3]
}
\)
consts = {
γ -> 
            1.4, 
(Subscript[R, g]) -> 
            8.314/32*^-3
};

\(\(setvars = {
P\_C == 
            3.5*^6, 
P\_1 == 
            6*^6, 
\(m\& . \) ==  .44, 
ω == 2  π\ 85000/60\  (*\ 
            85000\ RPM\ in\ rad/s\ *) , 
T\_1 == 
            200, 
r\_2 == 
            3.254*^-3, 

 (*\ 
            Guesses\ for\ now\ *) 
r\_3 == 
            3*^-3, 
β\_2 == 
            0
};\)
\)
\(\(fulleq = 
        Join[eqthermo /. constants, setvars];\)
\)
\(getVars[x_] := 
      Complement[
        Variables[\((x /. constants)\) /. {Equal -> List, 
              ArcSin -> List, Rule -> List, Cos -> List, 
              Sin -> List, Sec -> List, 
              Power -> List}], {}]
\)
solveEq[x_, params_] := 
      Solve[Join[x  /.  consts,  setvars], 
        Complement[getVars[x], 
          params]]

thermosolution = solveEq[eqthermo, {}]

W //. eqgeomsolved  /.  thermosolution