RSS    

   Задача Дирихле - (реферат)

Задача Дирихле - (реферат)

Дата добавления: март 2006г.

    Задача Дирихле
    1. ПОСТАНОВКА ЗАДАЧИ
    Решить численно задачу Дирихле для уравнения Лапласа :
    (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 .

    2. ОПИСАНИЕ МЕТОДА РЕШЕНИЯ ПОСТАВЛЕННОЙ ЗАДАЧИ

Поставленная задача решается численно с помощью программы, реализующей метод сеток , разработанный для численного решения задачи Дирихле для уравнений эллептического типа.

Программа написана на языке C++ , в среде Borland C++ версии 3. 1. Ниже описан алгоритм работы этой программы.

1. На первом шаге область D дискретизируется. Она заменяется на область Dh путем разбиения области D параллельными прямыми по следующему правилу: yi=y0 ? ih, xj=x0 ? ih , i, j=0, 1, 2…. PР. Разбиение производится до тех пор, пока текущая прямая не будет лежать целиком вне области 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]|
    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    for(j=0; j    FillF_();
    CountDif();
    }
    cout    MakeFile();
    getchar();
    } //MAIN END

int IsLineFit(float par, char Axis) // does the line belong to the defined area {

    switch(Axis) (par    else return 0;
    case 'x': if (par    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    printf("\n");
    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    printf("\n");
    } // end CreateNet
    void RegArrayX() // Regulation of arrays X & Y
    {
    int LastUnreg = 13 ;
    while (LastUnreg ! = 0) {
    for(i=0; 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    } printf("\n");
    }
    void RegArrayY()
    {
    int LastUnreg = 33 ;
    while (LastUnreg ! = 0) {
    for(i=0; i    if (arrayY[i]>arrayY[i+1]) { tmp=arrayY[i];
    arrayY[i]=arrayY[i+1];
    arrayY[i+1]=tmp; }}
    LastUnreg=LastUnreg-1; }
    for (i=0; i    printf("\n"); } // End of Regulation
    void CrMtrD(void) //Create general Matrix
    {
    for(i=0; i    for(j=0; j    arrayP[i][j]. BelongsToDh=0; }
    for(i=0; i    for(j=0; 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("\n");
    }
    int IsFit(point Par) //does point belong to area D?
    {
    if ((Par. xx=1. 99) && (Par. yy>=Par. xx)
    && (Par. yy    else return 0;
    }
    void CreateDh_(void) //Create area Dh_
    {
    for(i=0; i    for(j=0; j    if (IsFit(arrayP[i][j])) arrayP[i][j]. BelongsToDh_=1;
    cout     cout     }
    void FillF(void) // calc function F(x, y) at area Dh_
    {
    for(i=0; i    for(j=0; j    if (arrayP[i][j]. BelongsToDh_==1)
    arrayP[i][j]. F=arrayP[i][j]. xx*pow(arrayP[i][j]. yy, 2);
    else arrayP[i][j]. F=0;
    }
    int IsInner(int i, int j) //Is point inner?
    {
    if ((arrayP[i-1][j]. BelongsToDh_==1) &&
    (arrayP[i+1][j]. BelongsToDh_==1) &&
    (arrayP[i][j+1]. BelongsToDh_==1) &&
    (arrayP[i][j-1]. BelongsToDh_==1)) return 1;
    else return 0;
    }
    void CreateDh(void) //Create area Dh
    {
    for(i=0; i    for(j=0; j    if ((arrayP[i][j]. BelongsToDh_==1) &&
    IsInner(i, j))
    arrayP[i][j]. BelongsToDh=1;
    }
    void FillF_() //calc new appr. values of F
    {
    for(i=0; i    for(j=0; j    if (arrayP[i][j]. BelongsToDh==1)
    arrayP[i][j]. F_=(arrayP[i-1][j]. F+arrayP[i+1][j]. F+
    arrayP[i][j-1]. F+arrayP[i][j+1]. F)/4;
    else arrayP[i][j]. F_=0; }
    }
    void CountDif() // find maximal difference abs(F-F_)
    {
    k=0;
    for(i=0; i    for(j=0; j    { if (arrayP[i][j]. BelongsToDh==1) {
    diff[k]=fabs(arrayP[i][j]. F_-arrayP[i][j]. F);
    k++; }}
    E1=diff[0];
    for (k=1; k    if (diff[k]>E1) 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    for(j=0; j    f    " "    f. close() ;
    }
    4. ГРАФИКИ РЕШЕНИЙ
    РИС. 1 шаг h=0. 2
    РИС. 2 шаг h=0. 1
    5. ВЫВОД

Функция f(x, y) является неотрицательной в области D. Полученное решение лежит целиком над плоскостью XOY . Для данного решения выполняется принцип максимума.


Новости


Быстрый поиск

Группа вКонтакте: новости

Пока нет

Новости в Twitter и Facebook

                   

Новости

© 2010.