Текст подпрограммы и версий de59r_c.zip de59d_c.zip de60r_c.zip de60d_c.zip de57r_c.zip de57d_c.zip de58r_c.zip de58d_c.zip |
Тексты тестовых примеров tde59r_c.zip tde59d_c.zip tde60r_c.zip tde60d_c.zip tde57r_c.zip tde57d_c.zip tde58r_c.zip tde58d_c.zip |
Вычисление решения двухточечной краевой задачи для обыкноновенного линейного самосопряженного дифференциального уравнения второго порядка с разрывными коэффициентами на равномерной сетке с помощью однородной консервативной разностной схемы второго порядка точности.
Решается двухточечная краевая задача для линейного самосопряженного уравнения второго порядка
(1) d (k(x) dy/dx) / dx - g(x) y = - r(x) , xN < x < xK
с разрывными коэффициентами
(2) k(x) , g(x) , r(x) ∈ Q( 2 ) [xN, xK]
( Q ( 2 ) [xN, xK] - класс функций, кусочно - непрерывных на [xN, xK] вместе со вторыми производными ), удовлетворяющими неpавенствам
(3) k(x) ≥ c > 0 , g(x) ≥ 0 .
Hа концах отрезка [xN, xK] задаются произвольные комбинации следующих краевых условий
при x = xN (4) 1) y( xN ) = GN , (5) 2) k dy/dx - ( σN y - μN ) = 0 , σN ≥ 0 , при x = xK (6) 1) y( xK ) = GK , (7) 2) - k dy/dx - ( σK y - μK ) = 0 , σK ≥ 0 ;
при σN > 0 условие (5) (а при σK > 0 условие (7)) является условием третьего рода, при σN = 0 (или при σK = 0) - условием второго рода.
Решение вычисляется на равномерной сетке узлов, заданной на [xN, xK], с помощью однородной консервативной разностной схемы, имеющей на указанном классе коэффициентов (2) второй порядок точности, при этом предполагается, что точки разрыва коэффициентов k (x), g (x) и r (x) уравнения (1) совпадают с узлами этой сетки.
Самарский A.A., Теория разностных схем, "Hаука", M., 1977.
int de59r_c (real *ck, real *cg, real *cr, integer *nx, real *xn, real *xk, logical *gn1, real *gn, logical *gk1, real *gk, real *y, doublereal *alp, doublereal *bet, integer *ierr)
Параметры
ck - |
одномерный вещественный массив длины nx - 1,
содержащий на входе в подпрограмму значения
коэффициента k (x) уравнения, а именно:
i - й элемент ck (i) массива содержит значение
функции k (x) в точке x = xn + (i - 1/2) * (xk - xn)/(nx - 1); функция k (x) должна принимать положительные значения; |
cg, cr - |
одномерные вещественные массивы длины 2*nx - 2,
содержащие на входе в подпрограмму значения коэффициентов
g (x) и r (x), соответственно, а именно: cg (1) содержит предел справа функции g (x) в начале xn отрезка интегрирования, (2*i) - й и (2*i + 1) - й элементы cg (2*i), cg (2*i + 1) массива cg содержат, соответственно, пределы слева и справа функции g (x) в узле сетки x = xn + i*(xk - xn)/(nx - 1) (i = 1, 2, ..., nx - 2), (2*nx - 2) - й элемент cg (2*nx - 2) массива cg содержит предел слева функции g (x) в конце xk отрезка интегрирования, элементы массива cr содержат аналогичные значения, вычисленные для функции r (x); |
nx - | число узлов равномерной сетки, включая концы отрезка интегрирования, в которых требуется вычислить решение краевой задачи; nx > 1 (тип: целый); |
xn, xk - | концы отрезка интегрирования, в которых задаются граничные условия (4) или (5) и (6) или (7), соответственно; xk > xn (тип: вещественный); |
gn1 - | логическая переменная (или логическая константа), имеющая на входе в подпрограмму значение true, если в начале xn отрезка интегрирования задано граничное условие первого рода (4), и false в противном случае, т.е. когда в xn задано условие 2 - го или 3 - его рода (5); |
gn - |
если gn1 = true, то gn представляет значение
решения в начале xn отрезка интегрирования; если gn1 = false, то gn является именем одномерного массива длины 2, первый элемент gn (1) которого содержит значение σn, а второй элемент gn (2) - значение mn, входящие в граничное условие (5) (тип: вещественный); |
gk1 - | логическая переменная (или логическая константа), имеющая на входе в подпрограмму значение true, если в конце xk отрезка интегрирования задано граничное условие первого рода (6), и false в противном случае, т.е. когда в xk задано условие 2 - го (σk = 0) или 3 - го (σk > 0) рода (7); |
gk - |
если gk1 = true, то gk представляет значение
решения в конце отрезка интегрирования; если gk1 = false, то gk является именем одномерного массива длины 2, первый элемент gk (1) которого содержит значение σk, а второй элемент gk (2) - значение μk, входящие в граничное условие (7) (тип: вещественный); |
y - |
одномерный вещественный массив длины nx,
содержащий на выходе из подпрограммы решение
краевой задачи, вычисленное на равномерной
сетке; при этой i - й элемент y (i) этого
массива содержит приближенное значение решения в узле xn + (i - 1)*(xk - xn)/(nx - 1); |
alp - bet | одномерные вещественные рабочие массивы двойной точности длины nx - 1; |
ierr - | целая переменная, служащая для сообщения об ошибках, обнаруженных в процессе работы подпрограммы; при этом: |
ierr=65 - | если в граничном условии (5) коэффициент σn < 0; |
ierr=66 - | если в граничном условии (7) коэффициент σk < 0; |
ierr=67 - | если какое - нибудь значение коэффициента k (x), заданного в массиве ck, меньше или pавно 0; |
ierr=68 - | если какое - нибудь значение коэффициента g (x), заданного в массиве cg, меньше 0; |
ierr=69 - | если значение nx ≤ 1; |
ierr=70 - | если xk ≤ xn; |
во всех указанных случаях интегрирование уравнения прекращается. |
Версии
de60r_c - | вычисление решения двухточечной краевой задачи для обыкновенного линейного самосопряженного дифференциального уравнения второго порядка с разрывными коэффициентами на неравномерной сетке. |
de57r_c - | вычисление решения двухточечной краевой задачи для обыкновенного линейного самосопряженного дифференциального уравнения второго порядка с непрерывными коэффициентами на равномерной сетке. |
de58r_c - | вычисление решения двухточечной краевой задачи для обыкновенного линейного самосопряженного дифференциального уравнения второго порядка с непрерывными коэффициентами на неравномерной сетке. |
de57d_c - | то же, что и de57r_c, при этом все арифметические операции с вещественными числами выполняются с удвоенным числом значащих цифр. |
de58d_c - | то же, что и de58r_c, при этом все арифметические операции с вещественными числами выполняются с удвоенным числом значащих цифр. |
de59d_c - | то же, что и de59r_c, при этом все арифметические операции с вещественными числами выполняются с удвоенным числом значащих цифр. |
de60d_c - | то же, что и de60r_c, при этом все арифметические операции с вещественными числами выполняются с удвоенным числом значащих цифр. |
Первые операторы этих подпрограмм имеют следующий вид: int de60r_c (real *ck, real *cg, real *cr, integer *nx, real *x, logical *gn1, real *gn, logical *gk1, real *gk, real *y, doublereal *alp, doublereal *bet, integer *ierr) int de57r_c (R_fp fk, R_fp fg, R_fp fr, integer *nx, real *xn, real *xk, logical *gn1, real *gn, logical *gk1, real *gk, real *y, doublereal *alp, doublereal *bet, integer *ierr) int de58r_c(R_fp fk, R_fp fg, R_fp fr, integer *nx, real *x, logical *gn1, real *gn, logical *gk1, real *gk, real *y, doublereal *alp, doublereal *bet, integer *ierr) Параметр x в подпрограммах de60r_c и de58r_c - имя одномерного вещественного массива длины nx, в последовательных элементах которого, начиная с первого, размещаются в возрастающем порядке узлы неравномерной сетки, причем начало и конец отрезка интегрирования, в которых заданы граничные условия (4) - (7), хранятся, соответственно, в первом и последнем элементах x (1), x (nx) этого массива. Параметры fk, fg, fr в подпрограммах de57r_c и de58r_c - имена вещественных подпрограмм - функций float fk (float *x), float fg (float *x), float fr (float *x), вычисляющих значения коэффициентов уравнения k (x), g (x) и r (x), соответственно, при вещественном значении аргумента x = x. Остальные параметры указанных версий имеют такой же смысл, что и одноименные параметры подпрограммы de59r_c. При этом значение параметра ierr = 70 в подпрограммах de60r_c и de58r_c указывает на то, что узлы сетки расположены в массиве x не в возрастающем порядке, значения ierr = 67, 68 в п/п de58r_c и de57r_c - что значения коэффициентов k (x) и g (x) уравнения, вычисленные с помощью подпрограмм - функций fk и fg, не удовлетворяют требуемое неpавенство (3). B подпрограмме de60r_c в массиве ck размещаются значения коэффициента k (x) уравнения, вычисленные в "потоковых" точках, т.е. в середине между узлами сетки, а в массивах cg и cr - пределы слева и справа функций g (x) и r (x) в узлах неравномерной сетки. Kpоме того, каждая из подпрограмм de57r_c, de58r_c, de59r_c и de60r_c имеет версию de57d_c, de58d_c, de59d_c и de60d_c, соответственно, с теми же самыми эффектом и списком параметpов, выполняющую все операции над действительными числами с удвоенным числом значащих цифр. При этом параметры fk, fg, fr (в подпрограммах de57d_c,
de58d_c), являются именами подпрограмм - функций |
ck, cg - cr | (в подпрограммах de59d_c, de60d_c), |
xn, xk - | (в подпрограммах de57d_c, de59d_c), |
x - | (в подпрограммах de58d_c, de60d_c), |
gn, gk - y, alp bet | (во всех подпрограммах) |
должны иметь тип double. Тип остальных параметров не изменяется. В версиях подпрограмм de57r_c, de58r_c, de59r_c и de60r_c, формальные параметры alp и bet имеют тип double. |
Вызываемые подпрограммы
utde14_c - | подпрограмма выдачи диагностических сообщений при работе подпрограмм de57r_c, de58r_c, de59r_c и de60r_c. |
utde15_c - | подпрограмма выдачи диагностических сообщений при работе подпрограмм de57d_c, de58d_c, de59d_c и de60d_c. |
Замечания по использованию
Значение параметра nx, задаваемое при обращении к подпрограмме, должно быть достаточно большим, чтобы получить хорошее приближение численного решения к точному решению задачи. При интегрировании дифференциального уравнения (1) с непрерывными коэффициентами с помощью подпрограмм de57r_c, de58r_c, de57d_c, de58d_c значения коэффициента k (x) вычисляются в "потоковых" точках, т.е. в середине между узлами, а значения коэффициентов g (x) и r (x) - в узлах сетки последовательно слева направо (начало и конец xn, xk отрезка интегрирования также считаются узлами сетки). Поэтому эти продпрограммы можно использовать также и в случае, когда какой - нибудь из коэффициентов k (x), g (x) и r (x) (или, быть может, все) задан в виде массива значений в указанных точках. Эти значения должны быть расположены в последовательных элементах массива, начиная с первого. При этом: если в xn задано граничное условие 2 - го или 3 - го рода, то значения коэффициентов g (x) и r (x) должны быть заданы, начиная с первого узла xn сетки, а если 1 - го рода, то начиная со второго узла xn + (xk - xn)/(nx - 1); если в конце промежутка xk задано условие 1 - го рода, то значения коэффициентов g (x) и r (x) используются только до (nx - 1) - го узла включительно. B этом случае каждое обращение к соответствующей подпрограмме - функции вычисления коэффициента сводится к выборке значения очередного элемента массива. Значения параметров ck, cg, cr, nx, xk, x, gn1, gn, gk1, gk сохраняются. |
d (k(x) dy/dx) / dx = 0 , y(0) = 1 , y(1) = 0 , k(x) = 1 , x ≤ 0.55 , k(x) = 10 , x > 0.55 ; Точное решение: y = 1 - x/0.595 , x ≤ 0.55 , (1 - x) / 5.95 , x > 0.55 .int main(void) { /* Local variables */ extern int de59r_c(float *, float *, float *, int *, float *, float *, logical *, float *, logical *, float *, float *, double *, double *, int *); static int ierr, i__; static float y[21], cg[40], ck[20], gk, cr[40], gn, xk; static int nx; static float xn; static logical gk1, gn1; static double bet[20], alp[20]; nx = 21; xn = 0.f; xk = 1.f; gn1 = TRUE_; gn = 1.f; gk1 = TRUE_; gk = 0.f; for (i__ = 1; i__ <= 40; ++i__) { cg[i__ - 1] = 0.f; cr[i__ - 1] = 0.f; /* L200: */ } for (i__ = 1; i__ <= 10; ++i__) { /* L5: */ ck[i__ - 1] = 1.f; } for (i__ = 11; i__ <= 20; ++i__) { /* L15: */ ck[i__ - 1] = 10.f; } de59r_c(ck, cg, cr, &nx, &xn, &xk, &gn1, &gn, &gk1, &gk, y, alp, bet, &ierr); printf("\n %5i \n", ierr); printf("\n %16.7e \n", xn); printf("\n %16.7e \n", xk); printf("\n %5i \n", nx); for (i__ = 0; i__ <= 18; i__ += 3) { printf("\n %16.7e %16.7e %16.7e \n", y[i__], y[i__ + 1], y[i__ + 2]); } return 0; } /* main */ Результаты: ierr = 0 y(1) = 1.00000000000 + 00 y(11) = 9.09090909100 - 02 y(2) = 9.09090909089 - 01 y(12) = 8.18181818191 - 02 y(3) = 8.18181818179 - 01 y(13) = 7.27272727281 - 02 y(4) = 7.27272727271 - 01 y(14) = 6.36363636372 - 02 y(5) = 6.36363636362 - 01 y(15) = 5.45454545463 - 02 y(6) = 5.45454545453 - 01 y(16) = 4.54545454552 - 02 y(7) = 4.54545454545 - 01 y(17) = 3.63636363643 - 02 y(8) = 3.63636363636 - 01 y(18) = 2.72727272732 - 02 y(9) = 2.72727272727 - 01 y(19) = 1.81818181822 - 02 y(10) = 1.81818181818 - 01 y20) = 9.09090909109 - 03 Абсолютная погрешность во всех узлах сетки pавна 0.00000000000