|
Текст подпрограммы и версий de54r_c.zip , de54d_c.zip , de56r_c.zip , de56d_c.zip |
Тексты тестовых примеров tde54r_c.zip , tde54d_c.zip , tde56r_c.zip , tde56d_c.zip |
Вычисление решения двухточечной краевой задачи для линейного дифференциального уравнения второго порядка на равномерной сетке методом прогонки a.a.Абрамова с контролем точности.
Решается двухточечная краевая задача для линейного дифференциального уравнения второго порядка
(1) y '' + f(x) y ' + g(x) y = r(x)
c заданными на концах xN и xК отрезка интегрирования граничными условиями вида:
(2) aN y '(xN) + bN y(xN) = cN , a2N + b2N ≠ 0 (3) aK y '(xK) + bK y(xK) = cK , a2K + b2K ≠ 0
- методом прогонки А.А.Абрамова. Решение вычисляется на заданной на отрезке интегрирования равномерной сетке с заданной мерой погрешности EPS. Контроль точности по меpе погрешности заключается в следующем: если численное решение по абсолютной величине не меньше некоторой заданной константы P, то контроль точности ведется по относительной погрешности, иначе - по абсолютной.
А.А.Абрамов. Вариант метода прогонки. Ж. вычисл. матем. и матем. физ., 1961, 1, N 2, 349 - 351.
int de54r_c (real *ff, real *fg, real *fr, real *xn, real *xk,
real *cn, real *ck, integer *nx, real *x0, real *hx, logical *buldy,
real *hmin, real *eps, real *p, real *h, real *y, real *dy,
integer *nx1, integer *ierr)
Параметры
|
ff, fg - fr | имена вещественных подпрограмм - функций вычисления значений коэффициентов уравнения (1) f (x), g (x) и r (x), соответственно. Эти подпрограммы - функции должны иметь по одному паpаметpу, представляющему вещественное значение независимой переменной x; |
| xn, xk - | концы отрезка интегрирования, в которых заданы граничные условия (2) и (3), соответственно; xk > xn или xk < xn (тип: вещественный); |
| cn, ck - |
одномерные вещественные массивы длины 3, в
которых помещаются коэффициенты граничных
условий (2) и (3), т.е. cn (1) = an , cn (2) = bn , cn (3) = cn , ck (1) = ak , ck (2) = bk , ck (3) = ck ; коэффициенты граничных условий должны удовлетворять неpавенствам: a2n + b2n ≠ 0 и a2k + b2k ≠ 0 ; |
| nx - | число узлов равномерной сетки на отрезке интегрирования, в которых требуется вычислить решение задачи; nx ≥ 1 (тип: целый); |
| xo - | начальный узел равномерной сетки; xo должен находиться между xn и xk или совпадать с каким-нибудь концом (тип: вещественный); |
| hx - | шаг равномерной сетки; hx > 0, если xk > xn и hx < 0, если xk < xn (тип: вещественный); |
| buldy - | логическая переменная или константа, принимающая при входе в подпрограмму значение TRUE_, если вместе с решением задачи требуется вычислить и его производную, и FALSE_ - в противном случае, т.е. когда необходимо иметь только решение задачи; |
| hmin - | минимальное значение абсолютной величины шага интегрирования, котоpое разрешается использовать при численном решении воспомагательных задач Коши, к которым сводится данная краевая задача (тип: вещественный); |
| eps - | допустимая меpа погрешности, с которой тpебуется вычислить решение краевой задачи (тип: вещественный); |
| p - | граница перехода, используемая при оценке меры погрешности решения (тип: вещественный); |
| h - | вещественная переменная, содержащая начальное значение шага интегрирования, используемое при численном решении воспомагательных задач Коши; |
| nx1 - | целая переменная, указывающая на выходе из подпрограммы число узлов сетки, в которых вычислено решение данной краевой задачи; nx1 ≤ nx; |
| y, dy - | одномерные вещественные массивы длины не более nx1, в которых запоминаются значения pешения и его производной, вычисленные в nx1 узлах сетки, при этом элементы y (i) и dy (i) содержат значения решения и производной в узле xo + (i - 1) * h, причем производная вычисляется только тогда, когда buldy=true (таким образом, предполагается, что узлы сетки занумерованы в направлении от xn до xk); |
| ierr - | целая переменная, служащая для сообщения об ошибках, обнаруженных в процессе работы подпрограммы; при этом: |
| ierr=1 - | когда последний узел сетки xo + (n - 1) * h выходит за конец xk интервала интегрирования, а предпоследний узел xo + (n - 2) * h является внутренним узлом, т.е. лежит внутри интервала интегрирования; |
| ierr=3 - | когда несколько узлов сетки (больше одного) выходят за конец xk интервала интегрирования, а самый последний узел, принадлежащий отрезку интегрирования, является внутренним узлом; |
| ierr=2 - | когда один из узлов сетки совпадает с концом xk интервала интегрирования, а все последующие узлы не попадают в интервал интегрирования; |
| во всех указанных выше случаях решение и его производная вычисляются только в тех узлах сетки, которые принадлежат отрезку интегрирования, а в случае ierr = 1 и ierr = 3 они вычисляются также в конце xk отрезка и запоминаются в y (nx1) и dy (nx1), соответственно; |
| ierr=65 - | когда в граничном условии (2) на конце xn коэффициенты an = bn = 0; |
| ierr=66 - | когда в граничном условии (3) на конце xk коэффициенты ak = bk = 0; |
| ierr=67 - | когда начальный узел xo сетки не принадлежит отрезку интегрирования; |
| ierr=68 и | ierr=70 - когда решение краевой задачи не может быть вычислено с требуемой точностью eps при заданном минимальном значении hmin шага интегрирования, при этом ierr = 68 указывает, что точность не достигается при прямой прогонке Абрамова, а ierr=70 - при обратной прогонке; |
| при ierr = 65, 66, 67, 68, 70 интегрирование уравнения прекращается; в случае ierr = 68 и ierr = 70 интегрирование можно повторить обращением к подпрограмме с новыми значениями h и hmin; |
| ierr=69 - | когда данная краевая задача не может быть решена методом прогонки А.А.Абрамова. |
Версии
| de54d_c - | вычисление решения двухточечной краевой задачи для линейного дифференциального уравнения второго порядка на равномерной сетке методом прогонки А.А.Абрамова с повышенной точностью. При этом параметры ff, fg, fr, xn, xk, cn, ck, xo, hx, hmin, eps, p, h, y, dy и формальный параметp в подпрограммах - функциях ff, fg и fr должны иметь тип double. |
| de56r_c - | вычисление решения двухточечной краевой задачи для линейного дифференциального уравнения второго порядка на неравномерной сетке методом прогонки A.A.Абрамова с контролем точности. Первый оператор подпрограммы имеет вид: |
|
int de56r_c (real *ff, real *fg, real *fr, real *xn, real *xk, real *cn, real *ck, integer *nx, real *x, logical *buldy, real *hmin, real *eps, real *p, real *h, real *y, real *dy, integer *nx1, integer *ierr) Здесь: x - одномерный вещественный массив длины nx, представляющий узлы неравномерной сетки, в которых требуется вычислить решение задачи, при этом узлы сетки располагаются в x в направлении от xn до xk. Остальные параметры подпрограммы de56r_c имеют тот же смысл, что и одноименные параметры подпрограммы de54r_c. При этом выход из подпрограммы de56r_c со значением ierr = 67 указывает на то, что первый узел сетки x (1) не принадлежит интервалу интегрирования. |
| de56d_c - | вычисление решения двухточечной краевой задачи для линейного дифференциального уравнения второго порядка на неравномерной сетке методом прогонки a.a.Абрамова с повышенной точностью. Первый оператор подпрограммы имеет тот же вид, что и в подпрограмме de56r_c; при этом параметры ff, fg, fr, xn, xk, cn, ck, x, hmin, eps, p, h, y, dy и формальный параметр в подпрограммах - функциях ff, fg, fr должны иметь тип double. |
Вызываемые подпрограммы
| utde12_c - | подпрограмма выдачи диагностических сообщений при работе подпрограмм de54r_c и de56r_c. |
| utde13_c - | подпрограмма выдачи диагностических сообщений при работе подпрограмм de54d_c и de56d_c. |
| Kpоме того, подпрограммы de54r_c, de56r_c и подпрограммы de54d_c, de56d_c используют рабочие подпрограммы de54ru_c и de54du_c, соответственно. |
Замечания по использованию
|
B общем случае заданная точность приближенного решения не гарантируется. Значения параметров xn, xk, cn, ck, nx, xo, hx, buldy, hmin, eps, p, x сохраняются. |
y '' + x2 y ' + xy = 10e- x*x (2x2 - 1) (2 - X) , - 1≤ x≤ 1
y ' - 2y = 0 при x = - 1 ,
y ' + 2y = 0 при x = 1
( точное решение задачи y (x) = 10e- x*x )
int main(void)
{
/* Initialized data */
static float cn[3] = { 1.f,-2.f,0.f };
static float ck[3] = { 1.f,2.f,0.f };
/* Local variables */
extern int de54r_c(R_fp, R_fp, R_fp, float *, float *, float *, float *,
int *, float *, float *, logical *, float *, float *,
float *, float *, float *, float *, int *, int *);
static float hmin;
static int ierr;
extern float f_c(), g_c();
static float h__;
static int i__;
static float p;
extern float r_c();
static float y[11];
static logical buldy;
static float x0, dy[11], hx, xk, xn;
static int nx, nx1;
static float eps;
xn = -1.f;
xk = 1.f;
nx = 11;
x0 = xn;
hx = .01f;
p = 11.f;
hx = .010050251256281407f;
x0 = 1.f - hx * 10.f;
buldy = TRUE_;
hmin = 1e-10f;
eps = 1e-5f;
h__ = .01f;
de54r_c((R_fp)f_c, (R_fp)g_c, (R_fp)r_c, &xn, &xk, cn, ck, &nx, &x0, &hx,
&buldy, &hmin, &eps, &p, &h__, y, dy, &nx1, &ierr);
printf("\n %5i \n", ierr);
printf("\n %5i \n", nx1);
for (i__ = 1; i__ <= nx1; ++i__) {
printf("\n %16.7e \n", y[i__-1]);
}
for (i__ = 1; i__ <= nx1; ++i__) {
printf("\n %16.7e \n", dy[i__-1]);
}
return 0;
} /* main */
float f_c(float *x)
{
/* System generated locals */
float ret_val;
ret_val = *x * *x;
return ret_val;
} /* f_c */
float g_c(float *x)
{
/* System generated locals */
float ret_val;
ret_val = *x;
return ret_val;
} /* g_c */
float r_c(float *x)
{
/* System generated locals */
float ret_val;
/* Builtin functions */
double exp(double);
/* Local variables */
static float x2;
x2 = *x * *x;
ret_val = (float)exp((float)(-x2)) * 10.f * (x2 * 2.f - 1.f) * (2.f - *x);
return ret_val;
} /* r_c */
Результаты:
ierr = 0
nx1 = 11
y dy
4.45260617314 + 00 - 8.01021606015 + 00
4.37238323319 + 00 - 7.95378250048 + 00
4.29273839036 + 00 - 7.89518713214 + 00
4.21369299696 + 00 - 7.83450451556 + 00
4.13526765439 + 00 - 7.77180951660 + 00
4.05748221040 + 00 - 7.70717722605 + 00
3.98035575729 + 00 - 7.64068288024 + 00
3.90390663078 + 00 - 7.57240178309 + 00
3.82815240978 + 00 - 7.50240922954 + 00
3.75310991678 + 00 - 7.43078043029 + 00
3.67879521918 + 00 - 7.35759043838 + 00