ПОСТАНОВКА ЗАДАЧИ
Решить численно задачу Дирихле для уравнения Лапласа :
(x,y)??D ; u|Г=xy2=f(x,y) ;
область D ограничена линиями: x=2 , x=4 , y=x , y=x+4 ;
(x0, y0 ) = (3, 5) .
Следует решить задачу сначала с шагом по x и по y : h=0.2, потом с шагом
h=0.1 . Точность решения СЛАУ ?=0.01 .
ОПИСАНИЕ МЕТОДА РЕШЕНИЯ ПОСТАВЛЕННОЙ ЗАДАЧИ
Поставленная задача решается численно с помощью программы,
реализующей метод сеток , разработанный для численного решения задачи
Дирихле для уравнений эллептического типа.
Программа написана на языке C++ , в среде Borland C++ версии 3.1. Ниже
описан алгоритм работы этой программы.
1. На первом шаге область D дискретизируется. Она заменяется на область
Dh путем разбиения области D параллельными прямыми по следующему
правилу: yi=y0 ? ih, xj=x0 ? ih , i,j=0,1,2…. РР
Разбиение производится до тех пор, пока текущая прямая не будет лежать
целиком вне области D. Получается множество точек (xi,yj).
2. За область Dh принимают те точки множества (xi,yj) , которые
попали внутрь области D, а также дополняют это множество граничными
точками.
3.Во всех точках области Dh вычисляются значения функции
f(xi,yj) .
4. За область Dh* принимаются все внутренние точки области Dh, т.е.
удовлетворяющие требованию:
(xi,yj)?? Dh* , если (xi+1,yj) ? Dh , (xi-1,yj)?? Dh ,
(xi,yj+1)?? Dh , (xi,yj-1)?? Dh .
5. Во всех точках области Dh* вычисляется функция F(N)*[i,j] (
индекс N обозначает номер итерации, на которой производится вычисление):
F(N)*[i,j]=(f(xi+1,yj) + f(xi-1,yj) +
f(xi,yj+1)???f(xi,yj-1))/4
6. Теперь если max | F(N+1)*[i,j] – F(N)*[i,j]|? ,взятый по всем
точкам области Dh* ,то задача решена;
если нет , то выполнять шаг 5 ( пересчитывать функцию
F(N)*[i,j] через значения F(N-1)*[i,j]) до тех пор, пока не выполнится
указанное условие.
3.ТЕКСТ ПРОГРАММЫ
#include
#include
#include
#include
#include int i,j,k; // Variables float h,x,y,tmp,E1; struct point { float xx; float yy; int BelongsToDh_; int BelongsToDh; float F; float F_; } p0,arrayP[13][33]; float arrayX[13]; float arrayY[33]; float diff[500]; void CreateNet(void); // Procedure Prototypes int IsLineFit(float Param); void CrMtrD(void); void RegArrayX(); void RegArrayY(); void CreateDh_(); int IsFit(point Par); void FillF(); void CreateDh(); int IsInner(int i,int j); void FillF_(); void CountDif(); void MakeFile(); void main(void) //MAIN { clrscr(); p0.xx = 3; p0.yy = 5; h = 0.2; p0.BelongsToDh_=1; p0.BelongsToDh=1; CreateNet(); RegArrayX(); RegArrayY(); CrMtrD(); CreateDh_(); FillF(); CreateDh(); FillF_(); CountDif(); while (E1>=0.005) { for(i=0;i8.0) || (par4.0) return 1; else return 0; } } void CreateNet(void) // Creation of Net (area D ) { x = p0.xx; i=0; arrayX[i]=x; while (!IsLineFit(x,’x’)) { x += h; i++; arrayX[i] = x; } x = p0.xx-h; i++; arrayX[i]=x; while (!IsLineFit(x,’x’)) { x -= h; i++; arrayX[i] = x; } for (i=0;iarrayX[i+1]) {double tmp=arrayX[i]; arrayX[i]=arrayX[i+1]; arrayX[i+1]=tmp;}} LastUnreg=LastUnreg-1; } for (i=0;iarrayY[i+1]) { tmp=arrayY[i]; arrayY[i]=arrayY[i+1]; arrayY[i+1]=tmp;}} LastUnreg=LastUnreg-1; } for (i=0;i=1.99) && (Par.yy>=Par.xx) && (Par.yyE1) E1=diff[k];} } void MakeFile() { ofstream f; FILE *f1=fopen(“surf.dat”,”w1″); fclose(f1); f.open(“surf.dat”,ios::out,0); for(i=0;i
Нашли опечатку? Выделите и нажмите CTRL+Enter
