Историки об Елизавете Петровне: Елизавета попала между двумя встречными культурными течениями, воспитывалась среди новых европейских веяний и преданий...
Автоматическое растормаживание колес: Тормозные устройства колес предназначены для уменьшения длины пробега и улучшения маневрирования ВС при...
Топ:
Определение места расположения распределительного центра: Фирма реализует продукцию на рынках сбыта и имеет постоянных поставщиков в разных регионах. Увеличение объема продаж...
Организация стока поверхностных вод: Наибольшее количество влаги на земном шаре испаряется с поверхности морей и океанов...
Методика измерений сопротивления растеканию тока анодного заземления: Анодный заземлитель (анод) – проводник, погруженный в электролитическую среду (грунт, раствор электролита) и подключенный к положительному...
Интересное:
Принципы управления денежными потоками: одним из методов контроля за состоянием денежной наличности является...
Аура как энергетическое поле: многослойную ауру человека можно представить себе подобным...
Финансовый рынок и его значение в управлении денежными потоками на современном этапе: любому предприятию для расширения производства и увеличения прибыли нужны...
Дисциплины:
 2022-09-11 |
 26 |
из
5.00
|
Заказать работу |
|
|
|
|
Рассмотрим однородную эллиптическую задачу, определенную в области 
(рисунок 3) с краевыми условиями.
Электропроводность определяется следующим образом:
СЛАУ решаются с точностью 10-8.
Размеры включений одинаковы и составляют .
Количество включений (отдельно по каждой из координат и общее) и процентное соотношение общей площади всех включений ко всей области моделирования приведены в таблице 6.
Таблица 6 – Соответствие количества включений и занимаемой площади
| Количество включений | Занимаемая площадь, % |
| 4x4=16 | 4.00 |
| 5x5=25 | 6.25 |
| 6x6=36 | 9.00 |
| 7x7=49 | 12.25 |
| 8x8=64 | 16.00 |
| 9x9=81 | 20.25 |
| 10x10=100 | 25.00 |
| 11x11=121 | 30.25 |
| 12x12=144 | 36.00 |
| 13x13=169 | 42.25 |
| 14x14=196 | 49.00 |
Результаты моделирования приведены в таблице 7 и на рисунке 19.
Таблица 7 – Эффективное сопротивление при различных площадях включений и удельных электрических сопротивлений
| Занимаемая площадь, % | Удельное электрическое сопротивление | ||||||
| 1000 | 100 | 10 | 1 | 0.1 | 0.01 | 0.001 | |
| 4.00 | 1.076 | 1.074 | 1.063 | 1.000 | 0.923 | 0.902 | 0.900 |
| 6.25 | 1.121 | 1.119 | 1.100 | 1.000 | 0.883 | 0.852 | 0.835 |
| 9.00 | 1.180 | 1.176 | 1.147 | 1.000 | 0.836 | 0.794 | 0.789 |
| 12.25 | 1.253 | 1.249 | 1.206 | 1.000 | 0.784 | 0.731 | 0.725 |
| 16.00 | 1.344 | 1.337 | 1.278 | 1.000 | 0.728 | 0.662 | 0.655 |
| 20.25 | 1.474 | 1.464 | 1.376 | 1.000 | 0.675 | 0.605 | 0.597 |
| 25.00 | 1.624 | 1.611 | 1.489 | 1.000 | 0.616 | 0.539 | 0.530 |
| 30.25 | 1.816 | 1.797 | 1.629 | 1.000 | 0.557 | 0.473 | 0.463 |
Рисунок 21 – Эффективное сопротивление при различных площадях включений и удельных электрических сопротивлений
Из рисунка 21 видно, что результаты моделирования соответствуют физике процесса (при увеличении стехиометрических объемных долей в гетерогенном компакте эффективное сопротивление растет, если включения непроводящие, и уменьшается при проводящих включениях). Подобные результаты были получены и для других регулярных гетерогенных сред [33].
|
|
Экономический расчёт
Что-то тут пишем
Заключение
В работе было рассмотрено решение двумерной эллиптической краевой задачи с контрастными коэффициентами многомасштабным методом конечных элементов. Для решения задачи был реализован алгоритм многомасштабного метода конечных элементов [10]. Для проверки результатов задача так же решалась с помощью классического метода конечных элементов [28]. Найденное обоими методами распределение скалярного потенциала практически совпадают (с точностью до 10-5). Однако многомасштабный метод позволяет решать задачи с более мелкими включениями за счет независимого построения сетки в каждом макроэлементе. Это так же позволяет эффективно распараллеливать метод (получено почти линейное ускорение).
Реализованный алгоритм не предполагает периодичности решения или мелкомасштабных включений и может быть использован при моделировании как геометрически неоднородных, так и разномасштабных включений. Для удобства анализа результатов моделирования был реализован визуализатор, позволяющий по полученному решению получить вид распределения скалярного потенциала (рисунки 11, 12) и векторного поля плотности тока (рисунки 12, 15).
Работа была представлена на студенческой научной конференции ФПМИ, где получила второе место, и на новосибирской межвузовской студенческой конференции «Интелектуальный потенциал Сибири», а так же была опубликована [34].
Список литературы
1. Allaire, G. A multiscale finite element method for numerical homogenization / G. Allaire, R. Brizzi // SIAM MMS. – 2005. – 4. – С. 790-812.
2. Banks, H.T. Homogenization of Periodically Varying Coefficients in Electromagnetic Materials / H.T. Banks [и др.] // SAMSI: Technical Report. – 2005. – N2005-2. – С. 1-20.
3. Bensoussan, A. Asymptotic analysis for periodic structures / A. Bensoussan, J.L. Lions, G. Paranicolaou. – Noth-Holland: New York, 1978. – 721 с.
|
|
4. Bochev, P. A Mathematical Framework for Multiscale Science and Engineering: The Variational Multiscale Method and Interscale Transfer Operators / P. Bochev [и др.] // SAND REPORT. – 2004. – 2004-2871. – С. 1-20.
5. Bondeson, A. Computational Electromagnetics / A. Bondeson, T. Rylander, P. Ingelström. – B.: Springer, 2000. – 231 с.
6. Durmaz, S. A numerical study on the effective thermal conductivity of composite materials / S. Durmaz. – IZMIR, 2004. – 240 с.
7. E, W. Multiscale modeling and computation / W. E, B. Engquist // Notices Amer. Math. Soc.. – 2003. – V. 50: № 9. – С. 1062-1070.
8. E, W. Analysis of the heterogeneous multiscale method for ellpiptic homogenization problems / W. E, P. Ming, P. Zhang // J. Am. Math. Soc.. – 2003. – vol. 18. – С. 121-156.
9. E, W. The heterogeneous multiscale methods / W. E, B. Engquist // Comm. Math. Sci. – 2003. – vol. 1, no. 1. – С. 87-132.
10. Efendiev Y. Multiscale Finite Element Methods: Theory and Applications / Efendiev Y., Hou T. Y. – B.: Springer, 2009. – 241 с.
11. Efendiev, Y. Accurate multiscale finite element methods for two-phase flow simulations / Y. Efendiev, V. Ginting, T. Hou // Journal of Computational Physics. – 2006. – vol. 220, no.1. – С. 155-174.
12. Gilbert, A.C. A Comparison of Multiresolution and Classical One-dimensional Homogenoztion Schemes / A.C. Gilbert // Applied and Computational Hamonic Analysis. – 1998. – vol. 5, no. 1. – С. 1-35.
13. Guenneau, S. Homogenization of 3D finite chiral photonic crystals / S. Guenneau, F. Zolla // Phisica B. – 2007. – vol. 394. – С. 145-147.
14. Jiang, P.-X. Experimental research on convection heat transfer in sintered porous plate channels / Pei-Xue Jiang // International Journal of Heat and Mass Transfer. – 2004. – №47. – С. 2086-2096.
15. Hou. T. Y. A Multiscale Finite Element Method for Elliptic Problems in Composite Materials and Porous Media / Hou. T. Y., X.-H. Wu // Journal of computational physics. – 1997. – 134. – С. 169-189.
16. Hou, T. Convergence of a multiscale finite element method for elliptic problems with rapidly oscillating coefficients / T. Hou, X.-H. Wu, Z. Cai // Mathematics of Computation. – 1999. – vol. 68: no. 227. – С. 913-943.
17. Saad, Y. Iterative Methods for Sparse Linear Systems / Y. Saad. – 2nd edition. – B.: SIAM Society for Industrial & Applied Mathematics, 2003. – 477 с.
18. Sangalli, G. Capturing small scales in elliptic problems using a residual-free bubbles finite element method [Электронный ресурс] / G. Sangalli. – Режим доступа: http://www.siam.org/journals/mms/x-x/41140.html.
19. Weinan, E. Principles of Multiscale Modeling / E. Weinan. – New Jersey: Princeton University, 2011. – 488 с.
20. Бахвалов, Н.С. Осреднение процессов в периодических средах. / Н.С. Бахвалов, Г.П. Панасенко. – М.: Наука, 1984. – 352 с.
21. Бахвалов, Н.С. Осреднение дифференциальных уравнений с частными производными с быстро осциллирующими коэффициентами. / Н.С. Бахвалов // ДАН. – 1975. – Т.221: №3. – С. 516-519.
22. Дульнев, Г.Н. Теплопроводность смесей и композитных материалов: Справочная книга / Г.Н. Дульнев, Ю.П. Заричняк. – Л.: Энергия, 1974. – 264 с.
23. Жиков, В.В. Усреднение дифференциальных операторов / В.В. Жиков, С.М. Козлов, О.А. Олейник. – М.: Физматлит, 1993. – 464 с.
24. Жиков, В.В. Усреднение и перколяция / В.В. Жиков, С.М. Козлов // УМН. – 1988. – Т.43: №4. – С. 169-170.
25. Иткина, Н.Б. Использования многомасштабного метода для моделирования процесса диффузии гетерогенных средах / Н.Б. Иткина, Э.П. Шурина, Ю.И. Шокин. – Новосибирск: Изд-во НГТУ, 2009. – 10 с.
|
|
26. Копысов, С.П. Об одном методе определения эффективных характеристик композитов с помощью вейвлет-преобразования / С.П. Копысов, Ю.А. Сагдеева // Интеллектуальные системы в производстве. – 2007. – номер 1. – С. 49-62.
27. Ольшанский, М.А. Лекции и упражнения по многосеточным методам / М.А. Ольшанский. – М.: Физматлит, 2005. – 168 с.
28. Соловейчик, Ю.Г. Метод конечных элементов для решения скалярных и векторных задач: учеб. пособие / Ю.Г. Соловейчик, М.Э. Рояк, М.Г. Персова. – Новосибирск: Изд-во НГТУ, 2007. – 896 с.
29. Федоренко, Р.П. Введение в вычислительную физику / Р.П. Федоренко. – М.: Изд-во МФТИ, 1994. – 528 с.
30. Чишкала, В.А. Теплопроводность композиционных материалов на основе гидрида алюминия и терморасширенного графита. / В.А. Чишкала [и др.] // Вестник Харьковского университета. – 2008. – № 794. – С. 81-84.
31. Шокин, Ю.И. Многомасштабные методы: учеб. пособие / Ю.И. Шокин, Э.П. Шурина, Н.Б. Иткина. – Новосибирск: Изд-во НГТУ, 2010. – 68 с. – (Современные многосеточные методы; Часть I).
32. Шурина, Э.П. Методы решения СЛАУ большой рзмерности / Э.П. Шурина, М.Ю. Баландин. – Новосибирск: Изд-во НГТУ, 2000. – 70 с.
33. Эпов, М.И. Численная гомогенизация электрических характеристик сред с контрастными мелкомасштабными включениями / М.И. Эпов, Э.П. Шурина, М.К. Артемьев // Доклады академии наук. – 2011. – том 442: №1. – С. 1-3.
34. Кутищева, А.Ю. Решение эллиптической краевой задачи с контрастными коэффициентами многомасштабным методом конечных элементов: тезисы доклада / А.Ю. Кутищева // Современные проблемы технических наук. – 2012. – Часть 1. – С. 11.
Приложение А. ММКЭ
class MultiscaleFiniteElementMethods
{
private:
char direct[50];
std::vector<int> _use; //нужна для предотващени записи в один участок памяти
int TEST;
void input(double sig);
double func(double x, double y, int EL, int el);
double gauss (double a, double b, double c, double d, int EL, int el);
void enter_kraev4(int number, int El, int bf, double value, RarefiedMatrix A_fine[]);
void MFE_fine_elem(int current_elem, RarefiedMatrix A_fine[]);
double mult(std::vector<double>&x, std::vector<double>&y);
void local_matrix(int elem, int El, double koef[4][4], double G_loc[4][4], double b_loc[4]);
void MMFE_local_matrix(int elem, RarefiedMatrix A_fine[], int size);
|
|
void MMFE();
protected:
std::vector <grid> all_grids; //хранит все сетки. Первый элемент: coarse-grid. Остальные fine-grid в порядке их глобальной нумерации
std::vector <EnterTheBoundaryConditions> all_edges; //все краевые
int G_1[2][2], M_1[2][2]; //матрицы одномерных элементов
/*блок функций, отвечающих за генерацию портрета матрицы*/
struct spisok //для формировании портрета матрицы
{
int elem;
spisok *next;
};
void add(int min, int max, spisok *j);
void AssembiyPortrait(int el, RarefiedMatrix *Matrix, spisok *j);
void CreationPortrait(int el, RarefiedMatrix *Matrix);
void dimensional_matrix();
void global(int number_elem, int El, RarefiedMatrix *Matrix, double G_loc[4][4], double b_loc[4]);
void local_kraev3(int number, int El, RarefiedMatrix *Matrix);
void local_kraev2(int number, int El, RarefiedMatrix *Matrix);
void enter_kraev1(int number, int El, RarefiedMatrix *Matrix);
double X1(double x, double xmin, double xmax);
double X2(double x, double xmin, double xmax);
double Y1(double y, double ymin, double ymax);
double Y2(double y, double ymin, double ymax);
double dX(int i, double hx);
double dY(int i, double hy);
public:
bool _print;
bool _oscillating;
bool _parallel;
RarefiedMatrix A_coarse;//глобальная матрица coarse grid
double sigma[2];
MultiscaleFiniteElementMethods();
void PrintFullMatrix(char *file, int I, RarefiedMatrix A);
void PrintGl(char *file1, char *file2, RarefiedMatrix A, int I, int nX, int nY);
double PrintXY(double x, double y, std::vector<double> U, int El);
double PrintMsXY(double x, double y, std::vector<double> U);
std::vector <double> MultiscaleFiniteElementMethods::PrintJ(double x, double y, int El, int el);
double EffectiveResistance();
//запуск решателя.
//1ый параметр - деректория, где находится сетка (мб NULL)
//2ой - деректоирия, куда генерировать сетку (мб NULL)
//3ий - имя файла, где находятся данные для генерации сетки (мб NULL)
//sig - электропроводность
void calc(char *directory, char *generation, char *file_name, int tes, double sig);
~MultiscaleFiniteElementMethods();
};
#include "MMFE.h"
MultiscaleFiniteElementMethods::MultiscaleFiniteElementMethods(){}
void MultiscaleFiniteElementMethods::calc(char *directory, char *generation, char *file_name, int tes, double sig, int n_size)
{
A_coarse.Start("not symmetric");
if(directory!= NULL) strcpy(direct, directory);
if(generation!= NULL)
{
strcpy(direct, generation);
//запустить генератор
TheGeneratorOfGrid(file_name, direct, n_size);
}
TEST = tes;
input(sig);
MMFE();
}
/*ввод всех данных*/
void MultiscaleFiniteElementMethods::input(double sig)
{
char file_name[50];
char f1[50], f2[50], f3[50], f4[50], f5[50];
/*читаем файлы по сетке*/
/*сначала всю coarse-grid*/
sprintf(file_name,"test%d/coarse/parametrs.txt", TEST);
all_grids.push_back (grid(file_name, "2d"));
sprintf(file_name,"test%d/coarse/xy.txt", TEST);
all_grids[0].ReadFromFile_xy (file_name);
sprintf(file_name,"test%d/coarse/nvtr.txt", TEST);
all_grids[0].ReadFromFile_nvtr (file_name);
sprintf(file_name,"test%d/coarse/nvkat.txt", TEST);
all_grids[0].ReadFromFile_nvkat (file_name);
sprintf(f1,"test%d/coarse/kraev1.txt", TEST);
sprintf(f2,"test%d/coarse/kraev2.txt", TEST);
sprintf(f3,"test%d/coarse/kraev3.txt", TEST);
sprintf(f4,"test%d/coarse/kraev4.txt", TEST);
all_edges.push_back (EnterTheBoundaryConditions(f1,f2,f3,f4,f4, TEST));
/*читаем все по fine-grid*/
for(int i = 1; i <= all_grids[0].elem; i++)
{
sprintf (file_name,"test%d/fine/%d/parametrs.txt", TEST, i);
all_grids.push_back (grid(file_name, "2d"));
sprintf (file_name,"test%d/fine/%d/xy.txt", TEST, i);
|
|
all_grids[i].ReadFromFile_xy (file_name);
sprintf (file_name,"test%d/fine/%d/nvtr.txt", TEST, i);
all_grids[i].ReadFromFile_nvtr (file_name);
sprintf (file_name,"test%d/fine/%d/nvkat.txt", TEST, i);
all_grids[i].ReadFromFile_nvkat (file_name);
sprintf(f1,"test%d/fine/%d/kraev1.txt", TEST, i);
sprintf(f2,"test%d/fine/%d/kraev2.txt", TEST, i);
sprintf(f3,"test%d/fine/%d/kraev3.txt", TEST, i);
sprintf(f4,"test%d/fine/%d/kraev4.txt", TEST, i);
sprintf(f5,"test%d/fine/%d/kraev4_nvkat.txt", TEST, i);
all_edges.push_back (EnterTheBoundaryConditions(f1,f2,f3,f4,f5, TEST));
}
/*задаем значения коэфициэнтов*/
sigma[0] = sig;//включение
sigma[1] = 1.; //скелет
}
/*печатает матрицу в плотном виде*/
void MultiscaleFiniteElementMethods::PrintFullMatrix(char *file, int I, RarefiedMatrix A)
{
FILE *fmatr;
int temp,k;
fmatr = fopen(file, "w");
for(int i = 0; i < all_grids[I].node; i++)
{
for(int j = 0; j < all_grids[I].node; j++)
{
if(i == j)
{
if(A.diag[i] < 0) fprintf(fmatr, "%.2lf ", A.diag[i]);
else fprintf(fmatr, " %.2lf ", A.diag[i]);
}
if(i < j)
if(!A.getTypeMatrix())//матрица несимеетричная
{
temp=0;
for(k=A.ia[j]; k<A.ia[j+1] && temp==0; k++)
if(A.ja[k]==i)
{
temp=1;
if(A.au[k]<0)
fprintf(fmatr, "%.2lf ", A.au[k]);
else
fprintf(fmatr, " %.2lf ", A.au[k]);
}
if(temp==0) fprintf(fmatr, " xxxx ", 0.0);
}else{fprintf(fmatr, " xxxx ", 0.0);}
if(i > j)
{
temp=0;
for(k=A.ia[i]; k<A.ia[i+1] && temp==0; k++)
if(A.ja[k]==j)
{
temp=1;
if(A.al[k]<0)
fprintf(fmatr, "%.2lf ", A.al[k]);
else
fprintf(fmatr, " %.2lf ", A.al[k]);
}
if(temp==0) fprintf(fmatr, " xxxx ", 0.0);
}
if(j == all_grids[I].node - 1) fprintf(fmatr, "\n");
}
}
fprintf(fmatr, "\n b: ");
for(int i = 0; i < all_grids[I].node; i++)
fprintf(fmatr, "%.2lf ", A.b[i]);
fprintf(fmatr, "\n U: ");
for(int i = 0; i < all_grids[I].node; i++)
fprintf(fmatr, "%.4lf ", A.U[i]);
fprintf(fmatr, "\n U*: ");
double temp1=0, temp2=0;
for(int i = 0; i < all_grids[I].node; i++)
{
fprintf(fmatr, "%.4lf ", all_edges[I].u(all_grids[I].xy[0][i],all_grids[I].xy[1][i]));
temp1+=(A.U[i]-all_edges[I].u(all_grids[I].xy[0][i],all_grids[I].xy[1][i]))
*(A.U[i]-all_edges[I].u(all_grids[I].xy[0][i],all_grids[I].xy[1][i]));
temp2+=A.U[i]*A.U[i];
}
fprintf(fmatr, "\n delta: %e ", sqrt(temp1/temp2));
fclose(fmatr);
}
/*составление файлов для рисунка поля*/
void MultiscaleFiniteElementMethods::PrintGl(char *file1, char *file2, RarefiedMatrix A, int I, int nX, int nY)
{
FILE *fp;
double maxU=A.U[0], minU=A.U[0];
fp = fopen(file1, "w");
fprintf(fp, "%d %d\n", all_grids[I].node, all_grids[I].elem);
fprintf(fp, "%d %d\n", all_grids[I].numX, all_grids[I].numY); //колличество КЭ по каждому направлению
for(int i=0; i < all_grids[I].node; i++)
{
fprintf(fp, "%.15lf\n", A.U[i]);
if(maxU < A.U[i]) maxU = A.U[i];
if(minU > A.U[i]) minU = A.U[i];
}
fclose(fp);
fp = fopen(file2, "w");
fprintf(fp, "%lf\n", 0.007); //размерность одного пикселя
fprintf(fp, "%.15lf\n", maxU); //максимальное значение (для нормировки)
fprintf(fp, "%.15lf\n", minU); //минимальное
fclose(fp);
}
double MultiscaleFiniteElementMethods::X1(double x, double xmin, double xmax)
{
return (xmax - x) / (xmax - xmin);
}
double MultiscaleFiniteElementMethods::X2(double x, double xmin, double xmax)
{
return (x - xmin) / (xmax - xmin);
}
double MultiscaleFiniteElementMethods::Y1(double y, double ymin, double ymax)
{
return (ymax - y) / (ymax - ymin);
}
double MultiscaleFiniteElementMethods::Y2(double y, double ymin, double ymax)
{
return (y - ymin) / (ymax - ymin);
}
double MultiscaleFiniteElementMethods::dX(int i,double hx)
{
switch(i)
{
case 1: return -1./hx;
case 2: return 1./hx;
}
}
double MultiscaleFiniteElementMethods::dY(int i, double hy)
{
switch(i)
{
case 1: return -1./hy;
case 2: return 1./hy;
}
}
double MultiscaleFiniteElementMethods::PrintXY(double x, double y, std::vector<double> U, int El)
{
double Xmin, Xmax, Ymin, Ymax;
int elem, i, j;
for(i = 0; i < all_grids[El].X.size() - 1; i++)
{
Xmin = all_grids[El].X[i];
Xmax = all_grids[El].X[i+1];
if(x >= Xmin && x <= Xmax)
break;
}
for(j = 0; j < all_grids[El].Y.size() - 1; j++)
{
Ymin = all_grids[El].Y[j];
Ymax = all_grids[El].Y[j+1];
if(y >= Ymin && y <= Ymax)
break;
}
elem = i + j * (all_grids[0].X.size() - 1);
double rez = 0;
rez += U[all_grids[El].nvtr[0][elem]] * X1(x, Xmin, Xmax) * Y1(y, Ymin, Ymax);
rez += U[all_grids[El].nvtr[1][elem]] * X2(x, Xmin, Xmax) * Y1(y, Ymin, Ymax);
rez += U[all_grids[El].nvtr[2][elem]] * X1(x, Xmin, Xmax) * Y2(y, Ymin, Ymax);
rez += U[all_grids[El].nvtr[3][elem]] * X2(x, Xmin, Xmax) * Y2(y, Ymin, Ymax);
return rez;
}
double MultiscaleFiniteElementMethods::PrintMsXY(double x, double y, std::vector<double> U_global)
{
//определяем какому КЭ из coarse принадлежит точка
double Xmin, Xmax, Ymin, Ymax, Xmin_c, Xmax_c, Ymin_c, Ymax_c;
int elem_coarse, i, j;
for(i = 0; i < all_grids[0].X.size() - 1; i++)
{
Xmin_c = all_grids[0].X[i];
Xmax_c = all_grids[0].X[i+1];
if(x >= Xmin_c && x <= Xmax_c)
break;
}
for(j = 0; j < all_grids[0].Y.size() - 1; j++)
{
Ymin_c = all_grids[0].Y[j];
Ymax_c = all_grids[0].Y[j+1];
if(y >= Ymin_c && y <= Ymax_c)
break;
}
elem_coarse = i + j * (all_grids[0].X.size() - 1);
//опеределяем какому КЭ из fine
int elem_fine;
for(i = 0; i < all_grids[elem_coarse+1].X.size() - 1; i++)
{
Xmin = all_grids[elem_coarse+1].X[i];
Xmax = all_grids[elem_coarse+1].X[i+1];
if(x >= Xmin && x <= Xmax)
break;
}
for(j = 0; j < all_grids[elem_coarse+1].Y.size() - 1; j++)
{
Ymin = all_grids[elem_coarse+1].Y[j];
Ymax = all_grids[elem_coarse+1].Y[j+1];
if(y >= Ymin && y <= Ymax)
break;
}
elem_fine = i + j * (all_grids[elem_coarse+1].X.size() - 1);
//значение базисных функций
double rez_fine[] = {0.,0.,0.,0.};
for (i = 0; i < 4; i++)
{
rez_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[0][elem_fine]] * X1(x, Xmin, Xmax) * Y1(y, Ymin, Ymax);
rez_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[1][elem_fine]] * X2(x, Xmin, Xmax) * Y1(y, Ymin, Ymax);
rez_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[2][elem_fine]] * X1(x, Xmin, Xmax) * Y2(y, Ymin, Ymax);
rez_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[3][elem_fine]] * X2(x, Xmin, Xmax) * Y2(y, Ymin, Ymax);
}
//значение на coarse билинейном базисе
double rez = 0;
for(i = 0; i < 4; i++)
rez += U_global[all_grids[0].nvtr[i][elem_coarse]] * rez_fine[i];
return rez;
}
std::vector <double> MultiscaleFiniteElementMethods::PrintJ(double x, double y, int elem_coarse, int elem_fine)
{
std::vector<double> J(2);
double Xmin, Xmax, Ymin, Ymax;
Xmin = all_grids[elem_coarse+1].xy[0][all_grids[elem_coarse+1].nvtr[0][elem_fine]];
Xmax = all_grids[elem_coarse+1].xy[0][all_grids[elem_coarse+1].nvtr[1][elem_fine]];
Ymin = all_grids[elem_coarse+1].xy[1][all_grids[elem_coarse+1].nvtr[0][elem_fine]];
Ymax = all_grids[elem_coarse+1].xy[1][all_grids[elem_coarse+1].nvtr[2][elem_fine]];
//значение базисных функций
double Jx_fine[] = {0.,0.,0.,0.}, Jy_fine[] = {0.,0.,0.,0.};
for (int i = 0; i < 4; i++)
{
Jx_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[0][elem_fine]] * dX(1, Xmax - Xmin) * Y1(y, Ymin, Ymax);
Jx_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[1][elem_fine]] * dX(2, Xmax - Xmin) * Y1(y, Ymin, Ymax);
Jx_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[2][elem_fine]] * dX(1, Xmax - Xmin) * Y2(y, Ymin, Ymax);
Jx_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[3][elem_fine]] * dX(2, Xmax - Xmin) * Y2(y, Ymin, Ymax);
Jy_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[0][elem_fine]] * X1(x, Xmin, Xmax) * dY(1, abs(Ymin - Ymax));
Jy_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[1][elem_fine]] * X2(x, Xmin, Xmax) * dY(1, abs(Ymin - Ymax));
Jy_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[2][elem_fine]] * X1(x, Xmin, Xmax) * dY(2, abs(Ymin - Ymax));
Jy_fine[i] += all_grids[elem_coarse+1].U[i][all_grids[elem_coarse+1].nvtr[3][elem_fine]] * X2(x, Xmin, Xmax) * dY(2, abs(Ymin - Ymax));
}
//значение на coarse билинейном базисе
double sig = sigma[all_grids[elem_coarse+1].nvkat[elem_fine]];
for(int i = 0; i < 4; i++)
{
J[0] += sig * A_coarse.U[all_grids[0].nvtr[i][elem_coarse]] * Jx_fine[i];
J[1] += sig * A_coarse.U[all_grids[0].nvtr[i][elem_coarse]] * Jy_fine[i];
}
return J;
}
double MultiscaleFiniteElementMethods::EffectiveResistance()
{
double I = 0, S = 0, U;
for(int elem_coarse = 0; elem_coarse < all_grids[0].elem; elem_coarse++)
for(int elem_fine = 0; elem_fine < all_grids[elem_coarse+1].elem; elem_fine++)
I += gauss(all_grids[elem_coarse+1].xy[0][all_grids[elem_coarse+1].nvtr[0][elem_fine]],
all_grids[elem_coarse+1].xy[0][all_grids[elem_coarse+1].nvtr[1][elem_fine]],
all_grids[elem_coarse+1].xy[1][all_grids[elem_coarse+1].nvtr[0][elem_fine]],
all_grids[elem_coarse+1].xy[1][all_grids[elem_coarse+1].nvtr[2][elem_fine]],
elem_coarse,
elem_fine);
S = abs (all_grids[0].xy[0][all_grids[0].nvtr[0][0]] - all_grids[0].xy[0][all_grids[0].nvtr[1][all_grids[0].numX-1]]);
U = abs(all_edges[0].EnterFunc1xy(0,0, 1) - all_edges[0].EnterFunc1xy(0,0, 3));
return S * U / I;
}
double MultiscaleFiniteElementMethods::func(double x, double y, int EL, int el)
{
std::vector<double> J;
J = PrintJ(x,y,EL,el);
return sqrt (J[0]*J[0] + J[1]*J[1]);
}
// Integrate function f (x, y, z) on the [a,b]x[c,d]x[e,f];
double MultiscaleFiniteElementMethods::gauss (double a, double b, double c, double d, int EL, int el)
{
std::vector < double > xj (4), q (4);
/* Points */
xj [0] = sqrt ((3.0 - 2.0 * sqrt (6.0 / 5.0)) / 7.0);
xj [1] = sqrt ((3.0 + 2.0 * sqrt (6.0 / 5.0)) / 7.0);
xj [2] = -xj [1];
xj [3] = -xj [0];
/* Weights */
q [0] = (18 + sqrt (30.0)) / 36.0;
q [1] = (18 - sqrt (30.0)) / 36.0;
q [2] = q [1];
q [3] = q [0];
double s = 0, xi, eta, theta, tmp = 0;
for (int i = 0; i < 4; ++i)
{
xi = (a + b + xj [i] * (b - a)) / 2.0;
for (int j = 0; j < 4; ++j)
{
eta = (c + d + xj [j] * (d - c)) / 2.0;
s += q[i] * q[j] * func (xi, eta, EL,el);
}
}
return s * (d - c) * (b - a) / 4.0;
}
void MultiscaleFiniteElementMethods::add(int min, int max, spisok *jaa)
{
spisok *t, *q;
jaa[max].elem++;
for(t=&jaa[max]; t!=NULL;)
{
if(t->elem==min && &jaa[max]!=t)
{
t=NULL;
jaa[max].elem--;
}
else{
if(t->next==NULL)
{
t->next=new spisok;
t=t->next;
t->elem=min;
t->next=NULL;
t=NULL;
}else{
if(t->elem<min && t->next->elem>min)
{
q=t->next;
t->next=new spisok;
t->next->elem=min;
t->next->next=q;
t=NULL;
}else{
t=t->next;
}
}
}
}
}
void MultiscaleFiniteElementMethods::AssembiyPortrait(int el, RarefiedMatrix *Matrix, spisok *jaa)
{
int temp=0;
for(int i = 1; i < all_grids[el].node+1; i++)
temp += jaa[i-1].elem;
Matrix->CreateMatrix(all_grids[el].node, temp);
Matrix->ia[0] = 0;
spisok *t, *del;
for(int i = 1; i < all_grids[el].node + 1; i++)
Matrix->ia[i] = Matrix->ia[i-1] + jaa[i-1].elem;
for(int i = 0, j = 0; i < all_grids[el].node; i++)
{
Matrix->diag[i] = 0;
for(t = jaa[i].next; t!= NULL; j++)
{
Matrix->ja[j] = t->elem;
Matrix->al[j] = 0;
del = t;
t = t->next;
//удаляем всю структуру
delete del;
}
}
delete[] jaa;
}
void MultiscaleFiniteElementMethods::CreationPortrait(int el, RarefiedMatrix *Matrix)
{
int i, j;
spisok *jaa;
jaa = new spisok[all_grids[el].node];
for(i = 0; i < all_grids[el].node; i++)
{
jaa[i].next = NULL;
jaa[i].elem = 0;
}
for(i = 0; i < all_grids[el].elem; i++)
{
for(j = 0; j < 4; j++)
{
for(int l = 1; l + j < 4; l++)
{
if(all_grids[el].nvtr[j][i] > all_grids[el].nvtr[j+l][i])
add(all_grids[el].nvtr[l+j][i], all_grids[el].nvtr[j][i], jaa); //списки сразу упорядоченные
if(all_grids[el].nvtr[j][i] < all_grids[el].nvtr[l+j][i])
add(all_grids[el].nvtr[j][i], all_grids[el].nvtr[l+j][i], jaa);
}
}
}
AssembiyPortrait(el, &(*Matrix), jaa);
}
/*собственно реализация метода конечных элементов*/
void MultiscaleFiniteElementMethods::dimensional_matrix()
{
for(int i = 0; i < 2; i++)
for(int j = 0; j < 2; j++)
{
G_1[i][j] = abs(i-j) * (-2) + 1;
M_1[i][j] = abs(i-j) * (-1) + 2;
}
}
void MultiscaleFiniteElementMethods::local_matrix(int elem, int El, double koef[4][4], double G_loc[4][4], double b_loc[4])
{
int i,j;
double hx, hy, sig, temp;
double G_temp[4][4], M_loc[4][4];
dimensional_matrix();
hx=abs(all_grids[El].xy[0][all_grids[El].nvtr[0][elem]]-all_grids[El].xy[0][all_grids[El].nvtr[1][elem]]);
hy=abs(all_grids[El].xy[1][all_grids[El].nvtr[0][elem]]-all_grids[El].xy[1][all_grids[El].nvtr[2][elem]]);
if(El == 0)
{ //среднее значение на элементе
temp = 0.00006103515625 * (64 / all_grids[0].numX) * (64 / all_grids[0].numX);
sig = sigma[0] * temp + sigma[1] * (1-temp);
} else
sig = sigma[all_grids[El].nvkat[elem]];
G_temp[0][0] = sig * (G_1[0][0]*M_1[0][0]*hy/(hx*6.0) + G_1[0][0]*M_1[0][0]*hx/(hy*6.0));
G_temp[0][1] = G_temp[1][0] = sig * (G_1[0][1]*M_1[0][0]*hy/(hx*6.0) + G_1[0][0]*M_1[0][1]*hx/(hy*6.0));
G_temp[0][2] = G_temp[2][0] = sig * (G_1[0][0]*M_1[0][1]*hy/(hx*6.0) + G_1[0][1]*M_1[0][0]*hx/(hy*6.0));
G_temp[0][3] = G_temp[3][0] = sig * (G_1[0][1]*M_1[0][1]*hy/(hx*6.0) + G_1[0][1]*M_1[0][1]*hx/(hy*6.0));
G_temp[1][1] = sig * (G_1[1][1]*M_1[0][0]*hy/(hx*6.0) + G_1[0][0]*M_1[1][1]*hx/(hy*6.0));
G_temp[1][2] = G_temp[2][1] = sig * (G_1[1][0]*M_1[0][1]*hy/(hx*6.0) + G_1[0][1]*M_1[1][0]*hx/(hy*6.0));
G_temp[1][3] = G_temp[3][1] = sig * (G_1[1][1]*M_1[0][1]*hy/(hx*6.0) + G_1[0][1]*M_1[1][1]*hx/(hy*6.0));
G_temp[2][2] = sig * (G_1[0][0]*M_1[1][1]*hy/(hx*6.0) + G_1[1][1]*M_1[0][0]*hx/(hy*6.0));
G_temp[2][3] = G_temp[3][2] = sig * (G_1[0][1]*M_1[1][1]*hy/(hx*6.0) + G_1[1][1]*M_1[0][1]*hx/(hy*6.0));
G_temp[3][3] = sig * (G_1[1][1]*M_1[1][1]*hy/(hx*6.0) + G_1[1][1]*M_1[1][1]*hx/(hy*6.0));
for(i = 0; i < 4; i++)
for(j = 0; j < 4; j++)
{
G_loc[i][j] = 0.0;
if(El == 0) //т.е.элемент из coarse-grid
for(int k = 0; k < 4; k++)
G_loc[i][j] += koef[i][k] * G_temp[k][j];
else G_loc[i][j] = G_temp[i][j];
}
M_loc[0][0] = M_1[0][0] * M_1[0][0] * hx * hy / 36.0;
M_loc[0][1] = M_loc[1][0] = M_1[0][1] * M_1[0][0] * hx * hy / 36.0;
M_loc[0][2] = M_loc[2][0] = M_1[0][0] * M_1[0][1] * hx * hy / 36.0;
M_loc[0][3] = M_loc[3][0] = M_1[0][1] * M_1[0][1] * hx * hy / 36.0;
M_loc[1][1] = M_1[1][1] * M_1[0][0] * hx * hy / 36.0;
M_loc[1][2] = M_loc[2][1] = M_1[0][1] * M_1[1][0] * hx * hy / 36.0;
M_loc[1][3] = M_loc[3][1] = M_1[1][1] * M_1[0][1] * hx * hy / 36.0;
M_loc[2][2] = M_1[1][1] * M_1[0][0] * hx * hy / 36.0;
M_loc[2][3] = M_loc[3][2] = M_1[0][1] * M_1[1][1] * hx * hy / 36.0;
M_loc[3][3] = M_1[1][1] * M_1[1][1] * hx * hy / 36.0;
for(i = 0; i < 4; i++)
{
b_loc[i] = 0;
if (El == 0)
for(j = 0; j < 4; j++)
//формируем локальную правую часть
b_loc[i] += M_loc[i][j] * all_edges[0].f(all_grids[El].xy[0][all_grids[El].nvtr[j][elem]],
all_grids[El].xy[1][all_grids[El].nvtr[j][elem]],
all_grids[El].nvkat[elem]);
}
}
void MultiscaleFiniteElementMethods::global(int number_elem, int El, RarefiedMatrix *Matrix, double G_loc[4][4], double b_loc[4])
{
int i,j,k,q;
//забисываем в глобальную матрицу локальные
for(i = 0; i < 4; i++)
{ //обрабатываю только нижний треугольник
for(j = 0; j < i; j++)
{
if (all_grids[El].nvtr[i][number_elem] > all_grids[El].nvtr[j][number_elem])
{
q=0;
for (k = Matrix->ia[all_grids[El].nvtr[i][number_elem]]; k < Matrix->ia[all_grids[El].nvtr[i][number_elem]+1] && q!= 1; k++)
if (Matrix->ja[k] == all_grids[El].nvtr[j][number_elem])
{
q=1;
Matrix->al[k] += G_loc[i][j];
if(!Matrix->getTypeMatrix())
Matrix->au[k] += G_loc[j][i];
}
}else{
q=0;
for (k = Matrix->ia[all_grids[El].nvtr[j][number_elem]]; k < Matrix->ia[all_grids[El].nvtr[j][number_elem]+1] && q!= 1; k++)
if (Matrix->ja[k] == all_grids[El].nvtr[i][number_elem])
{
q = 1;
Matrix->al[k] += G_loc[i][j];
if(!Matrix->getTypeMatrix())
Matrix->au[k] += G_loc[j][i];
}
}
}
Matrix->diag[all_grids[El].nvtr[i][number_elem]] += G_loc[i][i];
//запихиваем локальный вектор b_loc в глобальный b
Matrix->b[all_grids[El].nvtr[i][number_elem]] += b_loc[i];
}
}
/*учет краевых условий*/
//number - краевое по счету (как в файле задавали)
//el - номер КЭ в coarse grid
//bf - номер базисной функции (для 4х краевых)
void MultiscaleFiniteElementMethods::enter_kraev4(int number, int El, int bf, double value, RarefiedMatrix A_fine[])
{
double temp = 0., h;
#if _oscillating //осцилирующие границы
temp = value;
#else
if(all_edges[El].kraev4[1][number] == 0)
{
h = abs(all_grids[0].xy[1][all_grids[0].nvtr[0][El-1]]-all_grids[0].xy[1][all_grids[0].nvtr[2][El-1]]);
if(bf == 0) temp = (all_grids[0].xy[1][all_grids[0].nvtr[2][El-1]] - all_grids[El].xy[1][all_edges[El].kraev4[0][number]])/h;
if(bf == 2) temp = (all_grids[El].xy[1][all_edges[El].kraev4[0][number]] - all_grids[0].xy[1][all_grids[0].nvtr[0][El-1]])/h;
}
if(all_edges[El].kraev4[1][number] == 1)
{
h = abs(all_grids[0].xy[0][all_grids[0].nvtr[0][El-1]]-all_grids[0].xy[0][all_grids[0].nvtr[1][El-1]]);
if(bf == 2) temp = (all_grids[0].xy[0][all_grids[0].nvtr[1][El-1]] - all_grids[El].xy[0][all_edges[El].kraev4[0][number]])/h;
if(bf == 3) temp = (all_grids[El].xy[0][all_edges[El].kraev4[0][number]] - all_grids[0].xy[0][all_grids[0].nvtr[0][El-1]])/h;
}
if(all_edges[El].kraev4[1][number] == 2)
{
h = abs(all_grids[0].xy[1][all_grids[0].nvtr[0][El-1]]-all_grids[0].xy[1][all_grids[0].nvtr[2][El-1]]);
if(bf == 1) temp = (all_grids[0].xy[1][all_grids[0].nvtr[2][El-1]] - all_grids[El].xy[1][all_edges[El].kraev4[0][number]])/h;
if(bf == 3) temp = (all_grids[El].xy[1][all_edges[El].kraev4[0][number]] - all_grids[0].xy[1][all_grids[0].nvtr[0][El-1]])/h;
}
if(all_edges[El].kraev4[1][number] == 3)
{
h = abs(all_grids[0].xy[0][all_grids[0].nvtr[0][El-1]]-all_grids[0].xy[0][all_grids[0].nvtr[1][El-1]]);
if(bf == 0) temp = (all_grids[0].xy[0][all_grids[0].nvtr[1][El-1]] - all_grids[El].xy[0][all_edges[El].kraev4[0][number]])/h;
if(bf == 1) temp = (all_grids[El].xy[0][all_edges[El].kraev4[0][number]] - all_grids[0].xy[0][all_grids[0].nvtr[0][El-1]])/h;
}
#endif
//дальше кусок из 1ых краевых
int i, j;
//вместо A.b[edge.kraev1[number]] ставим значение краевого u_b[number]
A_fine[bf].b[all_edges[El].kraev4[0][number]] = temp;//u_q[edge.kraev1[1][number]];
//гауссом обнуляем столбец edge.kraev1[0][number] (фактически работаем с правой частью)
for (i = 0; i < all_grids[El].node; i++)
{
if (i == all_edges[El].kraev4[0][number])
{
A_fine[bf].diag[i]=1;
//мы должны обнулить строку edge.kraev1[number]
for (j = A_fine[bf].ia[i]; j < A_fine[bf].ia[i+1]; j++)
A_fine[bf].al[j] = 0;
}
else
{
if (i < all_edges[El].kraev4[0][number])//как бы учет верхнего треугольника
{
for (j = A_fine[bf].ia[all_edges[El].kraev4[0][number]]; j < A_fine[bf].ia[all_edges[El].kraev4[0][number]+1]; j++)
if (A_fine[bf].ja[j] == i)
A_fine[bf].b[i] -= A_fine[bf].al[j] * temp;
}
else
for (j = A_fine[bf].ia[i]; j < A_fine[bf].ia[i+1]; j++)
if (A_fine[bf].ja[j] == all_edges[El].kraev4[0][number])
{
A_fine[bf].b[i] -= A_fine[bf].al[j] * temp;
A_fine[bf].al[j] = 0;
}
}
}
}
void MultiscaleFiniteElementMethods::local_kraev3(int number, int El, RarefiedMatrix *Matrix)
{
int i,j, k, q;
double A_S3[2][2], b_S[2];
double h_x, temp_b;
h_x = abs(all_grids[El].xy[0][all_edges[El].kraev3[0][number]] - all_grids[El].xy[0][all_edges[El].kraev3[1][number]])
+ abs(all_grids[El].xy[1][all_edges[El].kraev3[0][number]] - all_grids[El].xy[1][all_edges[El].kraev3[1][number]]);
for(i = 0; i < 2; i++)
{
temp_b=0;
for(j = 0; j < 2; j++)
{
A_S3[i][j] = all_edges[El].betta3(number) * (h_x/6.0) * M_1[i][j]; temp_b += A_S3[i][j
}
b_S[i] = temp_b;
}
//после формировки локальных сразу переталкиваем их в глобальник
if(Matrix->getTypeMatrix())//симметричная
{
for(i = 0; i < 2; i++)
{ //обрабатываю только нижний треугольник
for(j = 0; j < i; j++)
{
if(all_edges[El].kraev3[i][number] > all_edges[El].kraev3[j][number])
{
q=0;
for(k = Matrix->ia[all_edges[El].kraev3[i][number]]; k < Matrix->ia[all_edges[El].kraev3[i][number]+1] && q!=1; k++)
if(Matrix->ja[k] == all_edges[El].kraev3[j][number])
{
q=1;
Matrix->al[k] += A_S3[i][j];
}
}else{
q=0;
for(k = Matrix->ia[all_edges[El].kraev3[j][number]]; k < Matrix->ia[all_edges[El].kraev3[j][number]+1] && q!=1; k++)
if(Matrix->ja[k] == all_edges[El].kraev3[i][number])
{
q=1;
Matrix->al[k] += A_S3[i][j];
}
}
}
Matrix->diag[all_edges[El].kraev3[i][number]] += A_S3[i][i];
//запихиваем локальный вектор b_S в глобальный A.b
Matrix->b[all_edges[El].kraev3[i][number]] += b_S[i];
}
}else{
for(i = 0; i < 2; i++)
{ //обрабатываю только нижний треугольник
for(j = 0; j < i; j++)
{
|
|
|
История развития хранилищ для нефти: Первые склады нефти появились в XVII веке. Они представляли собой землянные ямы-амбара глубиной 4…5 м...
Таксономические единицы (категории) растений: Каждая система классификации состоит из определённых соподчиненных друг другу...
Поперечные профили набережных и береговой полосы: На городских территориях берегоукрепление проектируют с учетом технических и экономических требований, но особое значение придают эстетическим...
Особенности сооружения опор в сложных условиях: Сооружение ВЛ в районах с суровыми климатическими и тяжелыми геологическими условиями...
© cyberpedia.su 2017-2024 - Не является автором материалов. Исключительное право сохранено за автором текста.
Если вы не хотите, чтобы данный материал был у нас на сайте, перейдите по ссылке: Нарушение авторских прав. Мы поможем в написании вашей работы!
, 