[matlab] Численное решение системы ДУ по разным переменным

[AHummel]

Добрый день! Обращаюсь к Вам за помощью в решении системы дифференциальных уравнений.
В процессе выполнения научной работы встала задача численного решения системы для расчёта интенсивности лазерного излучения. Сама по себе запись выглядит довольно просто (обозначил переменные и константы в привычном виде):

$$\frac{d\beta(t, z)}{dt} = C_{1}(1-\beta)K_{p} - C_{2}\beta K_{l} - C_{3}\beta$$

$$\frac{dK_{p}(t, z)}{dz} = C_{4}(1-\beta)K_{p}$$

$$\frac{dK_{l}(t, z)}{dz} = C_{5}\beta K_{l}$$

Проблема поставлена полностью: есть начальные условия для решения задачи Коши ( $K_{p}(t, z=0), K_{l}(t, z=0), \beta(z, t=0)$ ), пределы дифференцирования, число проходов пучка в кристалле, количество импульсов и пр. Беда в том, что первое уравнение временное, а второе и третье – пространственные. И стандартные методы вроде ode45 неприменимы. Я думал про пошаговое решение методом predict-correct, но не знаю, как поступать после первого шага, когда $\beta$ уже не равно нулю.
Буду безмерно благодарен за любые идеи по реализации.
Спасибо!