Задача Дирихле

численно задачу Дирихле для уравнения Лапласа : (x,y)D; u | Г =xy 2 =f(x,y) ; область D ограничена линиями: x=2 , x=4 , y=x , y=x+4 ; (x 0 , y 0 ) = (3, 5) . Следует решить задачу сначала с шагом по x и по y : h=0.2, потом с шагом h=0.1 . Точность решения СЛАУ =0.01 . 2.ОПИСАНИЕ МЕТОДА РЕШЕНИЯ ПОСТАВЛЕННОЙ ЗАДАЧИ Поставленная задача решается численно с помощью программы, реализующей

метод сеток , разработанный для численного решения задачи Дирихле для уравнений эллептического типа. Программа написана на языке C++ , в среде Borland C++ версии 1. Ниже описан алгоритм работы этой программы. 1. На первом шаге область D дискретизируется. Она заменяется на область D h путем разбиения области D параллельными прямыми по следующему правилу: y i =y 0  ih,

x j =x 0  ih , i,j=0,1,2…. PР. Разбиение производится до тех пор, пока текущая прямая не будет лежать целиком вне области D. Получается множество точек (x i ,y j ). 2. За область D h принимают те точки множества (x i ,y j ) , которые попали внутрь области D, а также дополняют это множество граничными точками. 3.Во всех точках области D h вычисляются значения функции f(x i ,y j ) .

4. За область D h * принимаются все внутренние точки области D h , т.е. удовлетворяющие требованию: (x i ,y j ) D h * , если (x i+1 ,y j )  D h , (x i-1 ,y j ) D h , (x i ,y j+1 ) D h , (x i ,y j-1 ) D h . 5. Во всех точках области D h * вычисляется функция

F (N) *[i,j] ( индекс N обозначает номер итерации, на которой производится вычисление): F (N) *[i,j]=(f(x i+1 ,y j ) + f(x i-1 ,y j ) + f(x i ,y j+1 )f(x i ,y j-1 ))/6. Теперь если max | F (N+1) *[i,j] - F (N) *[i,j]|< ,взятый по всем точкам области D h * ,то задача решена; если нет , то выполнять шаг 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;i<13;i++) for(j=0;j<33;j++) arrayP[i][j].F=arrayP[i][j].F_; FillF_(); CountDif(); } cout<<(0-arrayP[7][17].F_); MakeFile(); getchar(); } //MAIIsLineFit(float par,char

Axis) // does the line belong to the defined area { switch(Axis) { case 'y': if ((par>8.0) || (par<2.0)) return 1; else return 0; case 'x': if (par<1.9) return 1; else if (par>4.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;i<13;i++) { printf("%g ",arrayX[i]);

} printf(" "); y = p0.yy; i = 0; arrayY[i]=y; while (!IsLineFit(y y')) { y += h; i++; arrayY[i] = y; } y = p0.yy - h; i++; arrayY[i]=y; while (!IsLineFit(y y')) { y -= h; i++; arrayY[i] = y; } for(i=0;i<33;i++) { printf("%g ",arrayY[i]);} printf(" "); } // end CreateNet void RegArrayX() // Regulation of arrays X & Y { int LastUnreg = 13 ; while (LastUnreg != 0) { for(i=0;i<

LastUnreg-1;i++) { if (arrayX[i]>arrayX[i+1]) {double tmp=arrayX[i]; arrayX[i]=arrayX[i+1]; arrayX[i+1]=tmp;}} LastUnreg=LastUnreg-1; } for (i=0;i<13;i++) { printf("%g ",arrayX[i]); } printf(" "); } void RegArrayY() { int LastUnreg = 33 ; while (LastUnreg != 0) { 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<33;i++) { printf("%g ",arrayY[i]); } printf("

");} // End of Regulation void CrMtrD(void) //Create general Matrix { for(i=0;i<13;i++) for(j=0;j<33;j++) {arrayP[i][j].BelongsToDh_=0; arrayP[i][j].BelongsToDh=0;} for(i=0;i<13;i++) for(j=0;j<33;j++) { arrayP[i][j].xx=arrayX[i]; arrayP[i][j].yy=arrayY[j]; } // printf("%g %g",arrayP[12][0].xx,arrayP[12][0]. yy); // printf(" "); } int IsFit(point Par) //does point belong to area

D? { if ((Par.xx<=4) && (Par.xx>=1.99) && (Par.yy>=Par.xx) && (Par.yy<=Par.xx+4)) return 1; else return 0; } void CreateDh_(void) //Create area Dh_ { for(i=0;i<13;i++) for(j=0;j<33;j++)