Изучение структур данных и алгоритмов, используемых в МКЭ

МОРФ

НГТУ

Кафедра прикладной математики

 и информатики

Лабораторная Работа №1

Факультет: ПМИ

Группа:         ПМ-06

Выполнили:  Бейсенов К.В.

    Ким М.Ю.

Проверил:    Чернышев А.В.

Петров Р.В.

Новосибирск

2004

Цель работы:

Изучение структур данных и алгоритмов, используемых в МКЭ, изучение приёмов вычисления локальных матриц, изучение сборки глобальных матриц в разрежено-строчном формате и правой части СЛАУ.

Задание:

1. ;
2. Треугольные конечные элементы;
3. Базисные функции линейные;
4. Первые однородные краевые условия.

Исходные данные:

nvtr.dat	Номера вершин и номер материала для треугольников	n1 n2 n3 m
xy.dat	Координаты вершин треугольника	x y
inf2tr.dat	Исходные данные о количестве узлов, краевых условиях…
mu	Проницаемость материала	nмат знач : имя
toku	Значения тока	nмат знач : имя
ig.2d jg.2d	Портрет для глобальной матрицы

Выходные данные:

gg.dat di.dat	Значения элементов матрицы
pr.dat	Значения элементов вектора правой части

Ход работы:

Будем считать, что функция u на каждом треугольнике является линейной функцией координат x и y. Введем на треугольнике три линейные функции,, такие, что функция равна единице в вершине 1 и нулю в двух остальных вершинах, функция  равна единице в вершине 2 и нулю в двух остальных вершинах, а функция  равна единице в вершине3 и нулю в двух остальных вершинах.

Нетрудно убедиться, что для этих функций справедливы соотношения:

                                

Текст программы:

// Assembly.cpp : Defines the entry point for the console application.

//

#include "stdafx.h"

#include "stdlib.h"

#include "math.h"

#define MUo 4.0*3.14159265358979e-7                   

// Функция открытия и провнрки открытия файла

FILE *Fopen(char *name, char *mode)

{

FILE *file;

if( (file = fopen(name, mode)) == NULL )

{

printf( "The file '%s' was not opened\n", name);

exit(1);

}

return file;

}

// Функция закрытия и проверки закрытия файла

void Fclose(FILE *file, char *name)

{

if( fclose(file) )

printf( "The file '%s' was not closed\n", name);

}

// Построение локальной матрицы и вектора на заданном Конечном Элементе

void SetLocalMatrixVector(double LocalMatrix[][3], double LocalVector[3], int NumberFiniteElement,

int **tops_triangle, double **xy, double *conductivity, double *current)

{

double Alpha[3][3];

int n1 = tops_triangle[NumberFiniteElement][0] - 1;

int n2 = tops_triangle[NumberFiniteElement][1] - 1;

int n3 = tops_triangle[NumberFiniteElement][2] - 1;

int m  = tops_triangle[NumberFiniteElement][3];

double determinant;

Alpha[0][1] = xy[n2][1] - xy[n3][1]; Alpha[0][2] = xy[n3][0] - xy[n2][0];

Alpha[1][1] = xy[n3][1] - xy[n1][1]; Alpha[1][2] = xy[n1][0] - xy[n3][0];

Alpha[2][1] = xy[n1][1] - xy[n2][1]; Alpha[2][2] = xy[n2][0] - xy[n1][0];

determinant = fabs( Alpha[2][2] * Alpha[1][1] - Alpha[1][2] * Alpha[2][1] );

for(int i = 0; i < 3; i++ )

{

LocalVector[i] = ( determinant / 6.0 ) * current[m];

for(int j = 0; j < 3; j++ )

LocalMatrix[i][j] = ( 1.0 / ( 2.0*conductivity[m]*determinant ) ) *

(Alpha[i][1]*Alpha[j][1] + Alpha[i][2]*Alpha[j][2]);

}

}

int _tmain(int argc, _TCHAR* argv[])

{

double *conductivity;               // проводимость мю

double *current;                     // токи

double **xy;                         // координаты вершин треугольников

double *gg;                          // массив элеметов глоб матматрицы

double *di;                          // массив с диагональными элементами

double *pr;                          // правая часть

int    *ig;                         // массив индексов для gg

int    *jg;                         // массив индексов для gg

int    **tops_triangle;             // вершины треугольников n1 n2 n3 m

int    amount_node;                 // количество узлов без Первых Краевых Условий

int    amount_triangle;             // количество треугольников

int    amount_element;              // количество элементов в матррице (для массива jg)

int    amount_marginal;             // количество первых краевых условий

double LocalMatrix[3][3];           // локальная матрица

double LocalVector[3];               // локальный вектор

FILE   *file;

char   symbol;

int    index;

int    index_i;

int    index_j;

// Считываем файл "mu"

file = Fopen("mu", "r");

conductivity = new double[10];

if( conductivity == NULL )

{

printf("Insufficient available memory for array 'conductivity'\n");

return 0;

}

while( !feof(file) )

{

fscanf(file, "%d", &index);

fscanf(file, "%lf", &conductivity[index]);

conductivity[index] *= MUo;  

fscanf(file, "%c", &symbol);

while( (symbol != '\n') && (!feof(file)) )

fscanf(file, "%c", &symbol);

}

Fclose(file,"mu");

// Считываем файл "toku"

file = Fopen("toku", "r");

current = new double[1000];

if( current == NULL )

{

printf("Insufficient available memory for array 'current'\n");