Привиет народ, мне нужна задача:"Решение краевой задачи для дифура 2го порядка"
Нужно ее запрограммить в Паскале... может у кого она есть или ссылка на нее, то пожалуйста дайте....
Полная версия: Диф.ур 2 го порядка
Пусть на отрезке [p,q] задано диф. уравн. 2 порядка.
r(x) d^2y/dx^2 + s(x)dy/dx + t(x)y = f(x)
и краевые условия
(d1 dy/dx +g1y) |(x=p) = _f1
(d2 dy/dx +g2y) |(x=q) = _f2
функции r(x), s(x), t(x), f(x) заданны.
d1,d2,g1,g2,_f1,_f2 - известные числа.
требуется определить y=y(x).
Для численного решения, введем на отрезке [p,q] равномерную сетку:
x0 = p, xn = q,
h= (p-q)/n, x(k+1) - xk = h, k = 0,1.... n
h-шаг сетки
Для апроксимации производных в диф. уравнении будем использовать разностные соотношения
(dy/dx)xk = (yk+1 - yk-1)/2h
(d^2y/dx^2)xk = (yk+1 - 2yk + yk-1) / (h^2)
Для апросимации производных в краевых условиях будем использовать соотнощшения:
(dy/dx)p = (y1-y0)/h
(dy/dx)q = (yn-yn-1)/h
Окончательно, получим соедующие разностыне уравнения:
c0y0 - b0y1 = f0
-akyk-1 + ckyk - bkyk+1 = fk
-anyn-1 + cnyn = fn
здесь
c0 = g1 + d1/h, b0 = - d1/h
ak = S(xk)/2h - r(xk)/h^2, ck = -2r(xk)/h^2 +t(xk)
bk = - r(xk)/h^2 - s(xk)/2h, k = 1,...n-1
an=d2/h
cn=g2+d2/h
Записав систему в матричном виде, мы получим
Ax=f
где A будет трехдиагональная матрица...
и если (вроде обязательно), для нее будет выполненн критерий регулярности Адомара (то есть будут все параметры строчного преобладания положительны) то ее можно факторизовать
LU разложением, поулчим
LUx=f
теперь Ux=y
Ly=f это простая систама, и
Ux=y тоже...
на практике для 3-ч диаг, матрицы существует очень быстрый метод, в отличии от метода Гаусса - метод прогонки
вот программа на фортране, при желании можно транслировать в паскаль.
Subroutine progon (y,altha, beta, n, p, q)
dimension y(1), altha(1), beta(1)
h=(q-p)/n
!h - шаг сетки
x=p
call coef(x,p,q,h,a,b,c,f)
! опр. коэфф. в уравнении
altha(1)=b/c
beta(1) = f/c
do 5 i=1,n
x=x+h
call coef (x,p,q,h,a,b,c,f)
beta(I+1) = (f+a*beta(i))/(c-a*altha(i))
if (i.eq.n) go to 5
altha(i+1) = b/(c-a*altha(i))
5 continue
y(n+1) = beta(n+1)
!теперь прогоночные коэфф. определены
do 10 i=n,1,-1
y(i)=altha(i)*y(i+1)+beta(i)
10 continue
!вычислили знач. искомой функции
return
end
Subroutine coef (x,p,q,h,a,b,c,f)
call gran(d1,g1,fi1,d2,g2,fi2)
!gran определяет коэфф. в краевых условиях
if (x.eq.p) then
a=0.0
b=-d1/h
c=g1+d1/h
f=fi1
return
end if
if (x.gt.q-h/2) then
a=d2/h
b=0.0
c=g2+d2/h
f=fi2
return
end if
call func(x,r,s,t,f)
! func для точки Х опр. коэфф. R,S,T,F в исходной диф. уравнении
a=s/(2*h)-r/(h**2)
b=-r/(h**2) - s/(2*h)
c=-2*r/(h**2)+t
return
end
Subroutine gran(d1,g1,fi1,d2,g2,fi2)
d1=0.0
g1=1.
fi1=0.0
d2=0.0
g2=1
fi2=0.0
return
end
Subroutine Func (x,r,s,t,f)
r=1.0
s=0.0
t=-2.0
f=-3.0*sin(x)
return
end
program TEST
! .......
call progon (y,altha, beta, n, p, q)
end program test
в ней, P и Q - границы отрезка... N-число узловых точек.
массив Y опр. значения искомой функции в узловых точках, массивы alpha и beta опр. прогоночные коэф..
В программе gran и func задаются исходные данные
в примере заданно уравнение
d^2/dx^2 - 2y = -3sinx
с граничными условиями Дирихле.
уфф..
r(x) d^2y/dx^2 + s(x)dy/dx + t(x)y = f(x)
и краевые условия
(d1 dy/dx +g1y) |(x=p) = _f1
(d2 dy/dx +g2y) |(x=q) = _f2
функции r(x), s(x), t(x), f(x) заданны.
d1,d2,g1,g2,_f1,_f2 - известные числа.
требуется определить y=y(x).
Для численного решения, введем на отрезке [p,q] равномерную сетку:
x0 = p, xn = q,
h= (p-q)/n, x(k+1) - xk = h, k = 0,1.... n
h-шаг сетки
Для апроксимации производных в диф. уравнении будем использовать разностные соотношения
(dy/dx)xk = (yk+1 - yk-1)/2h
(d^2y/dx^2)xk = (yk+1 - 2yk + yk-1) / (h^2)
Для апросимации производных в краевых условиях будем использовать соотнощшения:
(dy/dx)p = (y1-y0)/h
(dy/dx)q = (yn-yn-1)/h
Окончательно, получим соедующие разностыне уравнения:
c0y0 - b0y1 = f0
-akyk-1 + ckyk - bkyk+1 = fk
-anyn-1 + cnyn = fn
здесь
c0 = g1 + d1/h, b0 = - d1/h
ak = S(xk)/2h - r(xk)/h^2, ck = -2r(xk)/h^2 +t(xk)
bk = - r(xk)/h^2 - s(xk)/2h, k = 1,...n-1
an=d2/h
cn=g2+d2/h
Записав систему в матричном виде, мы получим
Ax=f
где A будет трехдиагональная матрица...
и если (вроде обязательно), для нее будет выполненн критерий регулярности Адомара (то есть будут все параметры строчного преобладания положительны) то ее можно факторизовать
LU разложением, поулчим
LUx=f
теперь Ux=y
Ly=f это простая систама, и
Ux=y тоже...
на практике для 3-ч диаг, матрицы существует очень быстрый метод, в отличии от метода Гаусса - метод прогонки
вот программа на фортране, при желании можно транслировать в паскаль.
Код
Subroutine progon (y,altha, beta, n, p, q)
dimension y(1), altha(1), beta(1)
h=(q-p)/n
!h - шаг сетки
x=p
call coef(x,p,q,h,a,b,c,f)
! опр. коэфф. в уравнении
altha(1)=b/c
beta(1) = f/c
do 5 i=1,n
x=x+h
call coef (x,p,q,h,a,b,c,f)
beta(I+1) = (f+a*beta(i))/(c-a*altha(i))
if (i.eq.n) go to 5
altha(i+1) = b/(c-a*altha(i))
5 continue
y(n+1) = beta(n+1)
!теперь прогоночные коэфф. определены
do 10 i=n,1,-1
y(i)=altha(i)*y(i+1)+beta(i)
10 continue
!вычислили знач. искомой функции
return
end
Subroutine coef (x,p,q,h,a,b,c,f)
call gran(d1,g1,fi1,d2,g2,fi2)
!gran определяет коэфф. в краевых условиях
if (x.eq.p) then
a=0.0
b=-d1/h
c=g1+d1/h
f=fi1
return
end if
if (x.gt.q-h/2) then
a=d2/h
b=0.0
c=g2+d2/h
f=fi2
return
end if
call func(x,r,s,t,f)
! func для точки Х опр. коэфф. R,S,T,F в исходной диф. уравнении
a=s/(2*h)-r/(h**2)
b=-r/(h**2) - s/(2*h)
c=-2*r/(h**2)+t
return
end
Subroutine gran(d1,g1,fi1,d2,g2,fi2)
d1=0.0
g1=1.
fi1=0.0
d2=0.0
g2=1
fi2=0.0
return
end
Subroutine Func (x,r,s,t,f)
r=1.0
s=0.0
t=-2.0
f=-3.0*sin(x)
return
end
program TEST
! .......
call progon (y,altha, beta, n, p, q)
end program test
в ней, P и Q - границы отрезка... N-число узловых точек.
массив Y опр. значения искомой функции в узловых точках, массивы alpha и beta опр. прогоночные коэф..
В программе gran и func задаются исходные данные
в примере заданно уравнение
d^2/dx^2 - 2y = -3sinx
с граничными условиями Дирихле.
уфф..
Это текстовая версия — только основной контент. Для просмотра полной версии этой страницы, пожалуйста, нажмите сюда.