Численное решение уравнения теплопроводности и волнового уравнения — различия между версиями

Версия 18:53, 8 марта 2015

Скачать программу: Equation v6-9-mini release.zip Текст программы на языке JavaScript (разработчики Цветков Денис, Кривцов Антон): <toggledisplay status=hide showtext="Показать↓" hidetext="Скрыть↑" linkstyle="font-size:default"> Файл "Equations.js" <syntaxhighlight lang="javascript" line start="1" enclose="div"> window.addEventListener("load", main_equations, false); function main_equations() {

   // Предварительные установки

   var a0 = 1;                             // масштаб расстояния
   var t0 = 1;                             // масштаб времени
   var m0 = 1;                             // масштаб массы
   var T0 = 1;                             // масштаб температуры

   var k0 = 2 * Math.PI / t0;              // масштаб частоты
   var C0 = m0 * k0 * k0;                  // масштаб жесткости

   var m = 1 * m0;                 	    // масса
   var C = 1 * C0;                 	    // жесткость

   var omega = Math.sqrt(C / m);
   var c = a0 * omega;

   var tw = CONFIG.N_Wave * a0 / c / 4;
   var beta =  Math.PI * c * c * Math.sqrt(tw);

   var dx = 1 * a0;                        // шаг сетки по оси x

   var spf = CONFIG.spf;                   // steps per frame - сколько расчетов проходит за каждый кадр отображения
   var dt = CONFIG.dt * t0;                // шаг интегрирования по времени

   var c3 = init_canvas(Equations_canvas); // инициализация canvas

   // слайдер
   slider_dt.max = dt * 2;
   slider_dt.min = dt / 50;
   slider_dt.step = dt / 100;
   slider_dt.value = dt;
   slider_dt.oninput = function() {dt = parseFloat(slider_dt.value);};

   var heat, wave;
   function start_new_system() {
       heat = new Heat_model(T0, dx, dt, beta);
       wave = new Wave_model(a0, k0, m, c, dx, dt);
   }
   function start_new_system_stair() {
       var func = function(M, M_max, N, val) {return IC_stair(M, M_max, N, val,
           parseFloat(stair_x1_number.value), parseFloat(stair_x2_number.value))};
       CONFIG.IC_Heat = func;
       CONFIG.IC_Wave = func;
       start_new_system()
   }
   function start_new_system_Gauss() {
       CONFIG.IC_Heat = IC_Gauss;
       CONFIG.IC_Wave = IC_Gauss;
       start_new_system()
   }

   button_restart_stair.onclick = start_new_system_stair;
   button_restart_Gauss.onclick = start_new_system_Gauss;

   start_new_system();
   control();

   function control() {
       heat.calculate_steps(spf, dt);
       wave.calculate_steps(spf, dt);
       draw_interface(c3);
       if (checkbox_graph_01.checked) draw_graph(c3, heat.get_D(), heat.get_div(), 1, "#ff0000", "u");
       if (checkbox_graph_02.checked) draw_graph(c3, wave.get_D(), wave.get_div(), 1, "#00ff00", "u");
       requestAnimationFrame(control);
   }

   function draw_graph(c, D, D_div, mult_y, color, val) {
       // график
       c.ctx.strokeStyle = color;
       c.ctx.beginPath();
       c.ctx.moveTo(0, c.h - D[0].u * (c.h - 20));
       for (var i = 1; i < D.length - 1; i++) {
           c.ctx.lineTo(i / (D.length - 1) * c.w, c.h -  D[i][val] / D_div * mult_y * (c.h - 20));
       }
       c.ctx.stroke();
   }

   function draw_interface(c) {
       c.ctx.clearRect(0, 0, c.w, c.h);

       // вертикальная линия
       c.ctx.strokeStyle = "#bbbbbb";
       c.ctx.beginPath();
       c.ctx.moveTo(c.w / 2, 0);
       c.ctx.lineTo(c.w / 2, c.h);
       c.ctx.stroke();
   }

   function init_canvas(canvas) {
       canvas.onselectstart = function () {return false;};     // запрет выделения canvas

       var canv_obj = {};
       canv_obj.ctx = canvas.getContext("2d");                 // на context происходит рисование
       canv_obj.w = canvas.width;                              // ширина окна в расчетных координатах
       canv_obj.h = canvas.height;                             // высота окна в расчетных координатах

       return canv_obj;
   }

}

function Heat_model(T0, dx, dt, X) {

   var N = CONFIG.N_Heat + 2;                                  // количество узлов по оси x + 2 для ГУ
   var T_start = 1 * T0;                                       // начальная температура

   var T = make_mass(N);
   var IC = CONFIG.IC_Heat;
   IC(T, T_start, N, "u");                                     // начальные условия
   BC_mirrow(T, N);                                            // граничные условия

   var div = 1 * T0;                                           // делитель для функции рисования
   this.get_div = function() {return div};
   var t = dt / 2;                                             // суммарное время с начала расчета

   this.get_D = function() {return T};
   this.calculate_steps = calc_implicit_Euler_method;

   var alpha = new Array(N + 2), beta = new Array(N + 2);
   function calc_implicit_Euler_method(spf, dt) {
       for (var s = 0; s < spf; s++) {
           var A = X * dt / (dx * dx + 2 * X * dt);
           var B = dx * dx / (dx * dx + 2 * X * dt);

           alpha[N - 2] = 0;
           beta[N - 2] = T[N - 2].u;
           for (var i = N - 3; i >= 1; i--) {
               alpha[i] = A / (1 - alpha[i + 1] * A);
               beta[i] = (B * T[i + 1].u + beta[i + 1] * A) / (1 - alpha[i + 1] * A);
           }

           T[1].u = (B * T[1].u + beta[1] * A) / (1 - alpha[1] * A - A);
           for (var i = 2; i < N - 1; i++) {
               T[i].u = (A * T[i - 1].u + B * T[i].u + beta[i] * A) / (1 - alpha[i] * A);
           }

           t += dt;
       }
   }

   return this;

} function Wave_model(a0, k0, m, c, dx, dt) {

   var N = CONFIG.N_Wave + 2;                                  // количество узлов по оси x + 2 для ГУ
   var U_start = a0 * k0 / (2 * Math.PI);

   var U = make_mass(N);
   var IC = CONFIG.IC_Wave;
   IC(U, U_start, N, "u");                                     // начальные условия
   IC_zero(U, U_start, N, "v");                                // начальные условия для скоростей
   BC_mirrow(U, N);                                            // граничные условия

   var div = 1 * a0 * k0 / (2 * Math.PI);                      // делитель для функции рисования
   this.get_div = function() {return div};
   var t = dt / 2;                                             // суммарное время с начала расчета

   this.get_D = function() {
       if (CONFIG.smooth_wave) return get_smoothed(U, N, "u", CONFIG.averaging_quality_Wave);
       else return U;
   };
   this.calculate_steps = function(spf, dt) {
       var koeff1 = c * c / (dx * dx) / m;
       for (var s = 0; s < spf; s++) {
           for (var i = 1; i < N - 1; i++) {
               U[i].v += (U[i + 1].u - 2 * U[i].u + U[i - 1].u) * koeff1 * dt;
           }
           for (var i = 1; i < N - 1; i++) {
               U[i].u += U[i].v * dt;
           }
           t += dt;
       }
   };
   return this;

}

function get_smoothed(D0, N, val, avl_quality) {

   var D = [];
   var avl = (N - 2) / avl_quality;
   var s = 1;
   var avl_half = avl / 2;
   for (var i = avl_half + 1; i < N - 1 ; i += avl) {
       var av = 0;
       for (var j = i - avl_half; j < i + avl_half; j++) {
           av += D0[j][val];
       }
       D[s] = {};
       D[s][val] = av / avl;
       s++;
   }
   // для первой точки добавляется объект, чтобы график располагался правильно по центру
   // но сама первая точка не добавляется - она не усреднена
   D[0] = {};
   return D;

} </source> Файл "Config.js" <syntaxhighlight lang="javascript" line start="1" enclose="div"> var CONFIG = {

    fps:       60                              // количество кадров в секунду
   ,spf:       200                             // сколько раз происходит вычисление между кадрами

   ,dt: 0.005

   // начальные условия
   //      IC_zero:                все значения - нули
   //      IC_stair:               ступенька посередине экрана
   //      IC_half_stair:          половина ступеньки
   //      IC_stair_random:        ступенька, составленная из случайных значений
   //      IC_half_stair_random:   половина ступеньки из случайных значений
   //      IC_Gauss:               распределение Гаусса

// ,IC_Heat: IC_Gauss // ,IC_Wave: IC_Gauss

   ,IC_Heat:   function(M, M_max, N, val) {return IC_stair(M, M_max, N, val, 0.33, 0.66)}
   ,IC_Wave:   function(M, M_max, N, val) {return IC_stair(M, M_max, N, val, 0.33, 0.66)}

   ,N_Heat:    1000                            // количество частиц по оси x
   ,N_Wave:    1000                            // количество частиц по оси x

   ,averaging_quality_Wave: 200                // на сколько отрезков усреднения поделить N (число должно быть делителем N)
   ,smooth_wave: true                         // сглаживать волновое уравнение

}; </source> Файл "IC_BC_LIB_v1.js" <syntaxhighlight lang="javascript" line start="1" enclose="div"> // Аргументы в функции IC и BC: // M - массив (сформированный с помощью функции make_mass данной библиотеки) // M_max - максимальное значение, которое должна задать функция (иногда данное значение передается, // но не испозьзуется, сделано для того, чтобы у всех фунций был одинаковый набор аргументов, // и можно было бы быстро менять используемую функцию (например, из файла настроек)) // N - Размер массива // val - ключ, по которому будет организован доступ к данным. В программе можно будет получить данные через него, // например, так: M[0].val // Чтобы лучше разобраться - просмотрите код какой либо программы, использующей эту библиотеку.

/**

* Создает специальный массив для последующего задания в нем начальных и граничных условий
* N - Длина массива (нужно учитывать, что первый и последний элементы будут использованы для граничных условий))
* */

function make_mass(N) {

   var M = [];
   for (var i = 1; i < N - 1; i++)
       M[i] = {};
   return M;

}

function IC_zero(M, M_max, N, val) { // Заполняет массив нулями

   for (var i = 1; i < N - 1; i++) {
       M[i][val] = 0;
   }

}

// 0 <= start, end <= 1 // start <= end // чтобы получить функцию вида IC_stair(M, M_max, N, val), используйте замыкание // например: var func = function(M, M_max, N, val) {return IC_stair(M, M_max, N, val, 0.33, 0.66)} function IC_stair(M, M_max, N, val, start, end) { // Задает "ступеньку" высоты M_max посередине экрана

   for (var i = 1; i < N - 1; i++) {
       if (i >= N * start && i < N * end) M[i][val] = M_max;
       else M[i][val] = 0;
   }

} function IC_half_stair(M, M_max, N, val) { // Задает мгновенный подъем от нуля до максимального значения

   for (var i = 1; i < N - 1; i++) {
       if (i <  N / 2) M[i][val] = 0;
       else M[i][val] = M_max;
   }

} function IC_stair_random(M, M_max, N, val) { // Задает "ступеньку" из случайных значений от -M_max до M_max

   for (var i = 1; i < N - 1; i++) {
       if (i >= N / 3 && i < 2 * N / 3) M[i][val] = (2 * Math.random() - 1) * M_max;
       else M[i][val] = 0;
   }

} function IC_half_stair_random(M, M_max, N, val) {// Правая половина экрана будет заполнена случайными значениями от -M_max до M_max

   for (var i = 1; i < N - 1; i++) {
       if (i <  N / 2) M[i][val] = 0;
       else M[i][val] = (2 * Math.random() - 1) * M_max;
   }

} function IC_Gauss(M, M_max, N, val) { // Задает распределение Гаусса.

   for (var i = 1; i < N - 1; i++) {
       var x = (i - 1) / (N - 2);
       M[i][val] = Math.exp(-Math.pow(x - 0.5, 2) * 50);
   }

}

function BC_mirrow(M, N) { // Зеркальные граничные условия

   M[0] = M[1];
   M[N - 1] = M[N - 2];

} function BC_periodic(M, N) { // Периодические граничные условия.

   M[0] = M[N - 2];
   M[N - 1] = M[1];

} </source> Файл "Equations.html" <syntaxhighlight lang="html" line start="1" enclose="div"> <!DOCTYPE html> <html> <head>

   <meta charset="UTF-8" />
   <title>Equations</title>
   <script src="IC_BC_LIB_v1.js"></script>
   <script src="Config.js"></script>
   <script src="Equations.js"></script>

</head> <body>

   <canvas id="Equations_canvas" width="1000" height="500" style="border:1px solid #000000;"></canvas>

   Скорость расчета: <input type="range" id="slider_dt" style="width: 150px;">


   <input type="checkbox" id="checkbox_graph_01" checked/>— Уравнение теплопроводности

   <input type="checkbox" id="checkbox_graph_02" checked/>— Волновое уравнение


   <input type="button" id="button_restart_Gauss" value="Рестарт Гаусс"/>


   <input type="button" id="button_restart_stair" value="Рестарт ступенька"/> (0 <= x1, x2 <= 1; x1 < x2)

   x1 = <input type="number" id="stair_x1_number" min=0 max=1 step=0.001 value="0.33"/>

   x2 = <input type="number" id="stair_x2_number" min=0 max=1 step=0.001 value="0.66"/>

</body> </html> </source> </toggledisplay>