用NDSolve函数,需给出边界条件：
v = 0.9;
a = NDSolve[{Theta'[
z] == -(v - (1 - v) (2. Cos[Theta[z]] + Cos[Theta[z]]^2))/R[z],
R'[z] == -1./Tan[Theta[z]], Theta[10] == 10, R[10] == 10}, {Theta,
R}, {z, -10, 10}]
Plot[Evaluate[{Theta[z], R[z]} /. a], {z, -10, 10}]
蓝色是Theta的曲线,红色是R的曲线：